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Industrial  and  Systems  Engineering 
(ABSTRACT) 

This  research  develops  the  availability  function  for  a  continuously 
demanded  component  which  is  maintained  by  an  age  replacement  preventive 
maintenance  policy.  The  availability  function,  A{t),  is  a  function  of  time  and  is 
defined  as  the  probability  that  the  component  functions  at  time  t.  The 
component  is  considered  to  have  two  states:  operating  and  failed.  In  this  policy, 
the  component  is  repaired  or  replaced  at  time  of  failure.  Otherwise,  if  the 
component  survives  T  time  units,  a  preventive  maintenance  service  is 
performed.  T  is  known  as  the  age  replacement  period  or  preventive 
maintenance  policy.  The  component  is  considered  to  be  as  good  as  new  after 
either  service  action  is  completed. 

A  renewal  theory  approach  is  used  to  develop  A{t).  Past  research  has 
concerned  infinite  time  horizons  letting  analysis  proceed  with  limiting  values. 

This  research  considers  component  economic  life  that  is  finite.  The  lifetime, 
failure  service  time  and  preventive  maintenance  service  time  probability 
distributions  are  unique  and  independent.  Laplace  transforms  are  used  to 
simplify  model  development.  The  age  replacement  period,  T,  is  treated  as  a 


parameter  during  model  development.  The  partial  Laplace  transform  is 
developed  to  deal  with  truncated  random  time  periods.  A  general  model  is 
developed  in  which  the  resulting  availability  function  is  dependent  on  both 
continuous  time  and  T.  An  exact  expression  for  the  Laplace  transform  of  A(t,  T) 
is  developed. 

Two  specific  cases  are  considered.  In  the  first  case,  the  lifetime,  repair 
and  preventive  maintenance  times  are  unique  exponential  distributions.  This 
case  is  used  to  validate  model  performance.  Tests  are  performed  for  f-»0,  f-»oo 
and  for  times  in  between  these  extremes.  Results  validate  model  performance. 
The  second  case  models  the  lifetime  as  a  Weibull  distribution  with  exponential 
failure  repair  and  preventive  maintenance  times.  Results  validate  model 
performance  in  this  case  also.  Exact  infinite  series  for  the  partial  and  normal 
Laplace  transform  of  the  Weibull  distribution  and  survivor  function  are  presented. 

Research  results  show  that  the  optimum  infinite  time  horizon  age 
replacement  period  does  not  maximize  average  availability  for  all  finite  values  of 
component  economic  life.  This  result  is  critical  in  lifecycle  maintenance 
planning. 
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CHAPTER  1 


INTRODUCTION 


1.1  Background 

Maintenance  costs  represent  a  major  portion  of  the  total  production 
costs  in  industrial  environments.  Surveys  have  shown  the  average  percentage 
of  cost  attributed  to  factory  maintenance  activities  ranges  from  about  fifteen  to 
forty  per  cent  with  an  average  around  twenty-five  per  cent  [1,2].  Obviously,  an 
effective  maintenance  planning  program  serves  to  reduce  these  costs  thus 
increasing  bottom-line  profits.  The  key  to  effective  maintenance  planning  is 
minimizing  unplanned  downtime.  For  example,  an  article  appearing  in  Service 
A/ews  [3]  and  Jacobson  and  Arora  [4]  report  that  computer  downtime  has  a 
huge  impact  on  the  bottom-line  profits  of  major  businesses.  A  survey  detailed 
in  [3]  shows  that  the  cost  to  American  corporations  due  to  computer  downtime 
is  about  four  billion  dollars  annually.  The  net  impact  of  computer  downtime  is 
estimated  at  $78,191  per  hour  and  $330,000  per  outage.  An  extreme  case  is 
Wall  Street,  where  one  minute  of  computer  downtime  means  losses  measured 
in  millions  of  dollars  [3,  4]. 
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The  implication  is  clear:  minimize  unplanned  downtime  through 
effective  planned  maintenance.  Planned  or  scheduled  maintenance  is  referred 
to  as  preventive  maintenance.  Most  preventive  maintenance  actions  require  a 
system  or  component  to  be  non-operational.  However,  the  advantage  is  that 
these  actions  may  be  scheduled  when  the  impact  on  production  and 
operational  requirements  is  minimized.  For  instance,  preventive  maintenance 
actions  might  be  scheduled  during  a  night  shift  or  on  a  week-end  when 
production  requirements  are  not  high.  The  objective  is  to  minimize  unplanned 
downtime  during  peak  times  of  production  and  operations.  As  stated  in  the 
examples  above,  unplanned  non-operational  periods  during  these  peak  times 
can  be  very  costly. 


1.2  Preventive  Maintenance 


Preventive  maintenance  policies  may  be  thought  of  as  the  replacement, 
servicing,  overhaul,  etc.  of  one  or  more  functioning  components  in  a  device  at 
selected  (scheduled)  points  in  time.  Thus,  the  preventive  maintenance  (PM) 
action  is  defined  as  required  for  the  specific  situation.  Assumptions 
concerning  the  preventive  maintenance  action  are  discussed  later.  Effective 
maintenance  planning  must  also  deal  with  the  unplanned  downtime  that  will 
occur.  As  noted,  these  emergency  or  failure  repairs,  replacements,  etc.  may 
be  very  costly.  The  objective  of  a  preventive  maintenance  policy  is  to  avoid 
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the  high  costs  of  failure  service  (unplanned  downtime)  by  performing 
preventive  maintenance  on  functioning  (non-failed)  components  and  thus 
accrue  the  usually  lower  costs  of  preventive  maintenance  service  (planned 
downtime). 

There  are  three  common  types  of  preventive  maintenance  policies. 
These  are  age  replacement,  block  replacement  and  opportunistic  preventive 
maintenance  policies.  An  age  replacement  policy  specifies  a  component  age, 
T,  at  which  a  component  undergoes  a  preventive  maintenance  action  (planned 
downtime).  If  the  component  fails  before  time  T  (unplanned  downtime),  the 
component  receives  a  failure  maintenance.  For  the  purposes  of  the  research 
described  here,  the  component  is  considered  as  good  as  new  after  either  a 
failure  repair  or  a  preventive  maintenance  action.  Thus,  the  clock  tracking 
component  age  is  reset  after  either  a  failure  repair  or  preventive  maintenance 
action. 

A  block  replacement  policy  specifies  a  defined  time  interval,  luock-  The 
component  is  preventively  maintained  at  integer  multiples  of  luock  regardless 
of  the  age  of  the  component.  If  the  component  fails  between  the  block 
replacement  times,  the  component  receives  a  failure  maintenance.  Thus  this 
policy  does  not  consider  when  the  last  component  failure  service  occurred.  A 
common  example  of  the  block  replacement  policy  is  the  periodic  replacement 
of  stadium  lights.  In  this  case  the  set-up  costs  tend  to  be  high  and  thus  it  is 
more  cost  effective  to  replace  lights  at  given  intervals  regardless  of  their  age. 
Berg  and  Epstein  [6]  provide  a  concise  description  and  comparison  of  age. 
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block,  and  failure  replacement  policies.  In  this  case,  a  failure  replacement 
policy  only  performs  failure  repair/replacement.  No  preventive  maintenance  is 
performed. 

The  opportunistic  replacement  policy  may  follow  any  one  of  an  infinitely 
large  number  of  preventive  maintenance  policies.  In  this  sense,  it  covers  all 
other  policies  that  do  not  fit  the  above  definitions  of  age  replacement  and  block 
replacement.  An  example  is  the  servicing  of  an  automobile.  During  a 
scheduled  oil  change,  an  inspection  reveals  that  the  brakes  are  nearing  the 
end  of  their  service  life.  The  brakes  are  replaced  even  though  they  may  have 
lasted  an  additional  amount  of  time.  The  brake  replacement  was 
"opportunistic". 


1.3  Preventive  Maintenance  Implementation  Concerns 

Administratively,  the  age  replacement  policy  is  more  difficult  to 
implement  than  the  block  replacement  policy.  Each  component  of  each 
subsystem  or  system  must  be  tracked  and  scheduled  for  preventive 
maintenance  actions.  Block  replacement  requires  less  administrative 
overhead  since  the  components  may  be  “grouped"  for  block  preventive 
maintenance.  Depending  on  the  complexity  of  an  opportunistic  policy,  the 
administrative  overhead  may  range  from  small  to  large.  For  large  complex 
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systems,  block  replacement  may  be  easier  to  implement.  However,  most 
companies  have  computer-based  maintenance  tracking  systems  that  allow 
more  complex  preventive  maintenance  policies  to  be  effectively  executed.  The 
preventive  maintenance  scheduling  policies  may  be  implemented  and 
automatically  tracked  by  a  computer-based  maintenance  action  system. 

In  reality,  opportunistic  maintenance  policies  will  intertwine  with  age 
replacement  or  block  replacement  policies.  Preventive  maintenance 
performed  according  to  an  age  replacement  or  block  replacement  policy  will  at 
times  reveal  other  problems  that  require  maintenance.  These  problems  may 
be  corrected  at  these  “opportunistic"  times. 

The  research  reported  in  this  dissertation  concerns  the  age  replacement 
preventive  maintenance  policy.  The  policy  is  examined  for  single  component 
systems  or  systems  that  may  be  modeled  as  a  single  component. 

Unscheduled  component  failure  repairs  and  scheduled  preventive 
maintenance  actions  are  considered  to  be  "perfect".  These  service  actions 
renew  the  component;  thus  the  component  is  as  good  as  new  at  the  conclusion 
of  a  service  action. 
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1 .4  Component  Availability 


The  effectiveness  of  repairable  components  and  systems  is 
appropriately  measured  by  a  parameter  known  as  availability.  Availability, 
A(t),  (a  function  of  time)  is  defined  as  the  probability  that  a  system  or 
component  is  operational  at  time  t.  The  limiting  availability,  A  =  lim  A{t), 

t—^OO 

may  be  thought  of  as  the  proportion  of  time  the  system  is  operational  in  the 
long  term.  Thus  in  the  long  term,  the  limiting  availability,  A,  is  [1,  7]; 


A  = 


_ E[uptime] _ 

E[uptime]+E[dovmtime] 


Other  measures  involving  A(t)  include  [7]: 


(1)  Average  availability,  Aa„g(T): 

Kvg{r)  =  ^  fo  A(f)  dt . 


(2)  Limiting  average  availability: 


(1.1) 


A:^,=  lim  Kvg{r)  =  lim  ^  Jq  A{t)  dt . 

T—^oo  r— >oo  ^ 
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Interval  average  availability  may  be  found  for  an  interval  [ti,  T2]  using 


-  Tl)  =  A(f)  dt  . 


In  the  computer  manufacturing  area,  many  manufacturers  are 
responding  to  customer  requirements  by  guaranteeing  availability  levels  for 
their  systems  [5],  Often  the  purchase  and/or  service  contracts  designate 
specific  time  intervals  where  availability  levels  are  evaluated.  If  the  contract 
availability  levels  are  not  met  at  these  designated  times,  the  manufacturer 
provides  compensation  (rebate,  extended  warranty,  etc.)  to  the  customer.  In 
these  cases,  using  limiting  availability  as  a  planning  tool  may  not  provide  the 
expected  results  from  either  the  manufacturer's  or  customer's  point  of  view.  A 
more  appropriate  measure  may  be  the  average  availability  for  the  time  period 
under  consideration.  Furthermore,  evaluation  of  an  optimal  age  replacement 
preventive  maintenance  policy,  T*.  using  the  limiting  case  where  time 
approaches  infinity  may  not  provide  the  expected  results.  The  key  to  the 
proper  evaluation  is  obtaining  the  function  A{t),  where  A{t)  provides  the 
probability  the  system  or  component  operates  as  a  function  of  time. 
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1.5  Research 


The  focus  of  this  research  effort  is  the  age  replacement  preventive 
maintenance  policy.  Previous  work  [7,12]  in  this  area  has  concentrated  on 
obtaining  the  optimal  replacement  period,  T*,  to  minimize  a  cost  rate  function. 
Most  models  assume  that  failure  maintenance  cost  is  ci  and  a  preventive 
maintenance  action  cost  is  c^.  Furthermore,  research  [7,12]  has  shown  that  for 
Cl  >  C2.  an  optimal  T*  exists  while  for  ci  <  C2,  T*  goes  to  infinity  (i.e.  perform  no 
preventive  maintenance).  These  results  are  obtained  by  considering  the  limit 
as  time,  t,  approaches  infinity.  Additionally,  if  ci  and  C2  are  taken  to  be  the 
expected  failure  repair  time  and  the  expected  preventive  maintenance  time, 
then  the  results  shown  above  also  maximize  the  limiting  availability. 

It  is  emphasized  that  in  both  cases,  a  limiting  result  is  developed  by 
taking  the  limit  as  t  approaches  infinity.  When  applied  to  actual  systems  and 
components,  this  implies  an  infinite  time  horizon  for  system  or  component 
operational  use.  Most  components  and  systems  have  finite  time  horizons  over 
which  they  are  operated  and  repaired.  At  the  end  of  the  system  or  component 
life,  the  system  is  often  replaced  by  newer  technology.  Many  examples  of  this 
exist  in  the  computer  and  telecommunications  industries.  Other  studies  have 
evaluated  the  case  where  failure  repair  times  and  preventive  maintenance 
times  are  stochastic  but  are  assumed  to  be  independent  and  identically 
distributed  random  variables. 

The  objective  of  this  research  is  to  evaluate  age  replacement  times,  T*, 
that  maximize  average  availability  for  finite  time  intervals.  This  maximization  of 
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the  average  availability  proceeds  unconstrained.  The  results  obtained  are 
compared  to  the  optimal  age  replacement  period  obtained  from  the  limiting 
case. 

The  first  task  of  this  research  is  the  development  of  a  stochastic  model 
to  allow  the  formulation  of  the  function  A{t)  for  an  age  replacement  preventive 
maintenance  policy.  The  stochastic  model  is  based  on  the  assumption  that  the 
component  lifetime,  failure  repair  times,  and  preventive  maintenance  times  are 
unique  random  variables  having  known  and  unique  distribution  functions. 

Even  though  past  research  has  recognized  that  failure  repair  and  preventive 
maintenance  times  may  be  distinct  random  variables  with  distinct  distribution 
functions,  the  limiting  case  analysis  allows  the  mean  value  of  these 
distributions  to  be  used  in  the  model  development.  Thus,  the  model  is 
developed  to  provide  an  optimal  replacement  period,  T^,  for  the  infinite  time 
horizon.  However,  the  model  development  methodology  in  the  limiting  case 
does  not  allow  the  development  of  the  availability  function,  A(t).  The 
availability  function,  A{t),  allows  both  the  time  interval  case  as  well  as  the 
limiting  cases  to  be  evaluated.  The  availability  function,  A{t),  is  developed 
analytically  and  is  numerically  estimated  for  specific  lifetime,  failure  repair  time 
and  preventive  maintenance  time  distributions. 

The  second  task  includes  the  evaluation  of  these  models  for  specific 
component  lifetime,  failure  repair  time  and  preventive  maintenance  time 
distributions.  First,  the  exponential  case  is  evaluated.  Unique  exponential 


9 


distributions  are  designated  for  the  failure  repair  time  and  the  preventive 
maintenance  time. 

The  second  case  to  be  treated  considers  aging,  i.e.,  the  component's 
lifetime  distribution  has  an  increasing  failure  rate  (IFR).  This  means  that  given 
the  component  survives  until  time  t,  the  probability  that  a  component  fails  in 
the  next  instant  of  time  increases  as  t,  the  conditional  survival  time  increases. 
One  widely  accepted  lifetime  distribution  exhibiting  the  IFR  characteristic  is  the 
Weibull  distribution  with  shape  parameter  greater  than  one  (/3  >  1).  Therefore, 
the  Weibull  failure  case  is  studied  in  combination  with  exponential  failure 
repair  and  preventive  maintenance  actions.  These  random  service  periods 
represent  the  component  or  system  downtime  (non-operational)  while  the 
component  or  system  lifetime  distribution  function  and  the  age  replacement 
period  govern  the  system  uptime  (operational). 

Specific  expressions  are  developed  for  A{t)  in  each  of  the  cases 
discussed  above.  For  each  case,  T.  the  age  replacement  policy  that 
maximizes  the  average  availability  over  the  interval  [0,  t]  is  found  using  a 
numeric  search.  Time  interval  values  considered  include  the  component's 
mean  lifetime  and  characteristic  life  in  the  case  of  Weibull  lifetime  distribution. 
Furthermore,  analysis  is  accomplished  at  small  values  of  time  as  well  as  large 
values  of  time  to  verify  availability  model  performance  near  t  =  0  and  the 
limiting  case,  t-*oo.  Whenever  possible  analytical  results  have  been  obtained 
before  resorting  to  specific  numerical  techniques  and  examples. 
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1.6  Summary 


This  research  develops  the  component  availability  function,  A(t),  for  the 
age  replacement  preventive  maintenance  policy  where  component  lifetime, 
failure  repair  and  preventive  maintenance  times  are  stochastic  and  have 
known  but  unique  distribution  functions.  The  model  development  assumes  that 
service  actions,  whether  a  failure  repair  or  preventive  maintenance,  provide 
complete  component  or  system  renewal.  Thus  the  component  is  as  good  as 
new  at  the  conclusion  of  a  service  action.  T,  the  age  replacement  period,  is 
treated  as  a  parameter  in  the  development  of  A{t).  We  find  the  age 
replacement  period,  T,  that  maximizes  average  availability  over  an  interval  [0,t] 
(unconstrained)  through  numerical  integration.  Intervals  examined  include  the 
component's  mean  lifetime  as  well  as  the  characteristic  life  in  the  case  of  the 
Weibull  lifetime  distribution.  The  characteristic  life  for  the  Weibull  lifetime 
distribution  is  defined  as  the  point  in  time,  t,  where  there  is  a  cumulative 
probability  of  component  failure  occurring  at  or  before  time  t,  equal  to  0.632 
[95]. 

Both  constant  failure  rate  (CFR)  (no  aging  effects)  and  increasing  failure 
rate  (IFR)  (aging)  component  lifetime  distributions  are  considered. 
Exponentially  distributed  failure  repair  times  and  preventive  maintenance  times 
are  examined  for  an  exponentially  distributed  lifetime  as  well  as  for  Weibull 
lifetimes.  Model  validation  is  accomplished  at  values  of  t  near  zero  and  for 
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large  values  of  t,  to  approximate  the  limiting  case.  Final  results  are  compared 
to  the  limiting  case. 
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CHAPTER  2 


LITERATURE  REVIEW 


2.1  Introduction 


Over  the  past  thirty-five  years  there  has  been  a  continuing  and  growing 
interest  in  the  use  of  preventive  maintenance  to  increase  system  productivity, 
product  quality,  and  cost  effectiveness.  This  is  not  surprising  considering  the 
vast  industrial  complex  engaged  in  producing  goods  and  delivering  services. 
The  vast  majority  of  these  systems  and  components  require  periodic 
maintenance  of  some  type  to  maintain  production  efficiency  as  well  as  product 
and  service  quality.  With  the  exception  of  electronic  components  and  systems, 
these  systems  and  components  age  with  use  and  time.  Their  performance 
may  be  termed  deteriorating  with  respect  to  the  probability  of  failure  over  time. 
When  this  characteristic  is  coupled  with  the  large  cost  of  unscheduled 
downtime,  it  is  logical  that  research  pertaining  to  preventive  maintenance 
policies  has  produced  many  models  for  systems  and  components  exhibiting 
deterioration. 
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The  purpose  of  this  review  is  to  give  an  overview  of  the  type  of 
preventive  maintenance  models  produced  by  past  research.  Specific  attention 
is  given  to  models  dealing  with  age  replacement  policies.  These  preventive 
maintenance  policies  repair  or  replace  a  component  or  system  at  failure  or 
perform  a  preventive  maintenance  action  at  a  pre-specified  time  interval,  T. 
This  time  interval  is  usually  considered  the  age  of  the  component  since  these 
models  consider  complete  component  or  system  renewal  at  the  conclusion  of  a 
failure  maintenance  or  preventive  maintenance  action.  Hence,  T  is  known  as 
an  age  replacement  policy  or  period.  Discussions  concerning  model 
development  and  analysis  techniques  such  as  renewal  theory  and  Laplace 
transforms  are  included.  The  review  includes  a  short  discussion  concerning 
probability  modeling  for  component  lifetimes  as  well  as  for  the  time  required  for 
repair/replacement  and  preventive  maintenance  action.  Finally,  a  summary  of 
the  presentation  concludes  this  chapter. 


2.2  Preventive  Maintenance  Models 


Four  major  articles  review  the  research  performed  concerning 
preventive  maintenance  models.  These  include  the  1965  paper  by  McCall  [8], 
the  1976  paper  by  Pierskalla  and  Voelker  [9],  the  1977  paper  by  Lie,  et  al  [10], 
and  the  1989  paper  by  Valdez-Flores  and  Feldman  [1 1  ].  These  papers 
provide  a  chronological  review  of  the  research  performed  concerning 
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preventive  maintenance.  The  content  of  each  of  these  papers  is  discussed 
below  in  chronological  order. 

In  addition  to  the  articles  named  above,  bibliographic  references  may  be 
found  in  Osaki  and  Nakagawa  [52],  Sherif  and  Smith  [53]  and  Sherif  [54], 


2.2.1  Preventive  Maintenance  Research  Prior  to  1965 


In  "Maintenance  Policies  For  Stochastically  Failing  Equipment:  A 
Survey"  [8],  McCall  divides  maintenance  models  into  two  distinct  categories. 
These  are  preventive  maintenance  and  preparedness  models.  In  preventive 
maintenance  models,  the  component,  system  or  equipment  is  subject  to 
stochastic  failure  and  the  state  of  the  component  is  always  known  with 
certainty.  In  what  McCall  terms  preparedness  models,  the  component  is 
subject  to  stochastic  failure,  however,  the  state  of  the  component  is  not  known 
with  certainty.  McCall's  discussion  is  divided  by  whether  the  component's 
lifetime  failure  distribution  is  known  or  not  known  with  certainty.  The  research 
performed  in  this  dissertation  pertains  to  age  replacement  preventive 
maintenance  policies  for  which  the  component's  lifetime  probability  density 
function  is  known.  Thus  discussion  of  past  research  in  this  time  period  is 
limited  to  this  area. 

The  periodic  preventive  maintenance  policy  repairs  or  replaces 
components  at  failure  or  preventively  maintains  them  at  age  T.  This  is  the  age 
replacement  preventive  maintenance  policy  defined  above.  The  periodic 
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preparedness  policy  inspects  and/or  repairs/replaces  equipment  at  age  T.  In 
the  simplest  situation,  both  models  consider  the  component  to  be  in  one  of  two 
states  -  either  an  operational  state  or  a  failed  state.  The  preventive 
maintenance  model  detects  component  failure  immediately  while  the 
preparedness  model  detects  component  failure  only  at  the  time  of  an 
inspection  or  repair/replacement. 

The  simplest  and  one  of  the  best  known  preventive  maintenance 
policies  is  the  strictly  periodic  or  age  replacement  policy.  Recall  that  in  this 
policy  the  component  is  repaired/replaced  upon  failure  or  when  the  component 
reaches  an  age  replacement  policy  of  T.  The  model  considers  complete 
component  renewal  at  the  completion  of  either  service  action.  Thus,  the 
component  age  T  is  the  same  as  the  time  since  the  last  service  action  was 
completed.  Barlow  and  Hunter  [12]  as  well  as  Barlow  and  Proschan  [13]  have 
shown  this  type  of  policy  is  only  worthwhile  to  implement  in  the  situation  where 
components  exhibit  the  effects  of  aging  (strictly  increasing  failure  rate,  see  [7]) 
and  the  cost  of  repair/replacement  due  to  component  failure  (unscheduled 
downtime)  is  greater  than  the  cost  of  a  preventive  maintenance  action.  These 
costs  may  also  be  represented  as  the  time  required  to  perform  the  component 
service  actions  described. 

In  the  case  of  cost  interpretation,  the  objective  is  to  minimize  cost  while 
in  the  case  of  time  interpretation,  the  objective  is  to  maximize  the  proportion  of 
time  the  component  is  operational.  This  is  also  known  as  limiting  availability. 
Both  analyses  are  equivalent  and  an  infinite  time  span  is  considered.  Even 
though  the  costs  (times)  for  failure  maintenance  and  preventive  maintenance 
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are  considered  to  be  random  variables  with  known  distributions,  the  infinite 
time  span  assumption  allows  the  expected  values  of  these  quantities  to  be 
used.  Thus,  one  may  derive  an  optimal  age  replacement  period  that  minimizes 
the  cost  rate  or  maximizes  the  limiting  component  availability  for  an  infinite 
time  horizon. 

However,  the  analysis  does  not  allow  for  the  derivation  of  the  cost  or 
availability  as  a  function  of  time.  Therefore,  the  calculation  of  the  cost  or 
average  availability  over  a  specified  time  interval  may  not  be  found  since  the 
cost  rate  or  availability  function  is  not  available  for  integration.  Other  research 
work  referenced  for  the  age  replacement  policy  includes  Barlow  and  Proschan 
[13],  Campbell  [14],  Dean  [15],  Drenick  [16],  Kamins  and  McCall  [17],  Klein  and 
Rosenberg  [18],  Lotka  [19],  Morse  [20]  and  Welker  [21  and  22].  Furthermore, 
both  Barlow  and  Proschan  [13]  and  Derman  [23]  have  shown  that  for  infinite 
time  horizons,  a  strictly  periodic  replacement  policy,  i.e.,  non-random,  is 
optimal  over  random  replacement  policies. 

Additionally,  the  McCall  article  discusses  two  preventive  maintenance 
policies  for  systems  or  equipment  made  up  of  several  parts.  An  opportunistic 
model  developed  by  Radner  [24]  considers  a  two  component  system  consisting 
of  component  1  and  component  2.  The  components  are  considered  to  fail 
stochastically  and  independently.  Component  1  has  a  strictly  increasing  failure 
rate  (deterioration)  while  the  component  2  has  a  constant  failure  rate. 
Component  2  does  not  exhibit  deterioration  and  thus  is  modeled  using  the 
exponential  probability  density  function  for  its  lifetime.  The  model  does 
recognize  economies  of  scale  in  maintenance.  This  means  that  it  is  less  costly 
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to  replace  both  components  at  one  time  than  it  is  to  replace  each  component 
separately.  Radner  [24]  termed  these  model  conditions  as  stochastic 
independence  and  economic  dependence.  He  shows  that  the  optimal  policy  is 
characterized  by  two  decision  parameters,  n  and  N.  If  the  age  of  the  system,  x, 
is  less  than  n,  (x  <  n),  then  component  1  is  replaced/repaired  only  if  it  fails.  If 
n  <  X  <  N,  then  replace  component  1  if  either  part  fails.  If  x  >  N,  then 
replace/repair  component  1  at  once. 

A  block  replacement  policy  is  defined  by  McCall  [8].  In  this  policy, 
system  components  are  preventively  repaired/replaced  at  regular  non-random 
intervals  regardless  of  their  age.  The  model  assumes  each  component  is 
repaired  upon  failure  as  well.  Comparisons  of  the  operating  characteristics 
between  block  and  age  replacement  policies  may  be  found  in  Barlow  and 
Proschan  [7  and  13],  Boivard  [25],  Cox  [26],  Flehinger  [27,  28  and  29],  Hunter 
and  Proschan  [30],  and  Arrow,  et  al  [31]. 

An  age  replacement  model  for  finite  time  horizon,  (0,  t),  has  been 
researched  and  analyzed.  This  policy  is  known  as  a  sequential  preventive 
maintenance  policy.  This  policy  repairs  or  replaces  the  component  upon 
failure  or  at  the  age  replacement  interval,  T.  However,  it  differs  from  the 
previous  infinite  time  horizon  case  in  that  at  each  event  (whether  failure  or 
preventive  maintenance),  the  optimal  replacement  period,  T,  is  recalculated  for 
the  remaining  time  period.  Thus,  this  replacement  period  minimizes  the 
expected  cost  over  the  remaining  period  (or  maximizes  the  component 
availability  depending  on  the  interpretation).  Barlow  and  Proschan  [13]  prove 
that  the  expected  cost  resulting  from  following  an  optimal  sequential  policy 
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during  the  interval  (0,  t)  is  always  less  than  or  equal  to  the  expected  cost  of  the 
corresponding  optimal  periodic  policy.  They  also  prove  that  the  optimal 
sequence  is  non-random  and  illustrate  a  method  for  calculating  the  optimal 
intervals. 

The  periodic  preparedness  model,  as  discussed  in  McCall  [8],  and 
under  the  conditions  where  the  component's  lifetime  distribution  is  known  may 
concern  either  an  optimal  replacement  period  or  an  optimal  inspection  period. 
For  the  case  where  an  optimal  replacement  period  is  found,  the  component's 
condition  is  only  known  at  the  time  of  replacement  and  the  prescribed  service 
action  is  always  carried  out  regardless  of  the  component's  condition. 
References  for  research  in  this  area  include  Jorgenson  and  McCall  [32], 
Jorgenson  [33],  Savage  [34],  Solomon  and  German  [35]  and  Veinott  and 
Wagner  [36].  In  the  other  case,  an  optimal  inspection  period  is  found.  The 
difference  in  this  case  is  that  the  component  is  inspected  for  proper  operation 
and  repaired/replaced  only  if  found  to  be  in  a  failed  state.  References  for  work 
in  this  area  include  [32],  [33],  Barlow,  et  al  [37],  Coleman,  et  al  [38],  and 
Kamins  [39].  Kamins  [39]  and  Coleman,  et  al  [38]  have  also  extended  the 
model  to  include  imperfect  inspection. 


2.2.2  Preventive  Maintenance  Research  Between  1965  and  1975 


Two  papers  provide  excellent  summaries  of  the  research  accomplished 
concerning  preventive  maintenance  during  this  time  period.  The  first,  "A 
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Survey  of  Maintenance  Models:  The  Control  and  Surveillance  of  Deteriorating 
Systems"  by  Pierskalla  and  John  A.  Voelker  [9]  surveys  the  research 
accomplished  since  the  1965  paper  by  McCall  [8].  Contrasting  this  survey  is 
the  paper  by  Lie,  et  al  [10],  titled  "Availability  of  Maintained  Systems;  A  State- 
of-the-Art  Survey".  This  paper  provides  a  systematic  classification  of  the 
literature  relevant  to  availability.  One  specific  area  is  the  effect  of  preventive 
maintenance  on  availability.  An  excellent  discussion  of  the  various  definitions 
of  availability  is  also  given  in  this  paper. 

Lie,  et  al  [10]  classify  availability  depending  on  the  time  interval 
involved.  These  are:  (1)  instantaneous  availability.  A{t),  (2)  average  uptime  or 
average  availability  for  a  time  interval  [0,  r],  designated  Aavgij),  and  (3)  steady 
state  availability,  A.  These  quantities  are  defined  as  follows; 

( 1 )  A{t)  =  Probability  the  component  is  operational  at  time,  t 

(2)  Aaygir)  =  }£A{t)dt 

(3)  A  =  lim  A(t)  . 

>00 

These  definitions  may  be  found  in  Barlow  and  Proschan  [7]  and  Sandler  [40]. 
Barlow  and  Proschan  [7]  also  define  the  limiting  average  availability,  ^oo-  This 
measure,  A^  ,  is  mathematically  defined  as 
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They  show  that  ^oo  is  equivalent  to  A  when  the  limit  exists.  Previously,  limiting 
availability  was  often  referred  to  the  limiting  efficiency  (see  Barlow  and  Hunter, 
1960  [12]). 

Lie,  et  al  [10]  claim  that  the  most  appropriate  measure  of  availability 
depends  on  a  component's  or  system's  mission  and  its  condition  of  use.  They 
provide  the  following  examples  of  application: 

The  steady-state  availability  may  be  a  satisfactory 
measure  for  systems  which  are  to  be  operated 
continuously,  for  example  a  detection  radar  system.  The 
average  uptime  [average  availability]  may  the  most 
satisfactory  measure  for  systems  whose  usage  is  defined 
by  a  duty  cycle,  for  example,  a  tracking  radar  system 
which  is  called  upon  only  after  an  object  has  been 
detected,  and  is  expected  to  track  continuously  during  a 
given  time  period.  For  systems  which  are  required  to 
perform  a  function  at  any  random  time,  the  instantaneous 
availability  may  be  the  most  satisfactory  measure.  An 
example  of  such  a  system  is  a  data  processing  system 
used  in  air  traffic  control  which  is  called  upon  to  process 
flight  paths,  and  then  remain  idle  for  a  length  of  time. 

Lie,  et  al  [10]  also  note  that  in  general,  several  types  of  preventive 
maintenance  policies  are  possible.  These  are  age  replacement,  block 
replacement,  random  periodic  replacement  and  sequential  determined 
replacement.  Numerous  references  are  cited  for  these  types  of  preventive 
maintenance.  They  note  that  the  earliest  research  on  planned  replacement 
was  by  Campbell  [14]. 
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Maintenance  models  for  this  period  are  extensively  surveyed  by 
Pierskalla  and  Voelker  [9].  This  paper  includes  only  those  models  which 
involve  an  optimal  decision  to  procure,  inspect,  repair  and/or  replace  a  single 
component  exhibiting  deterioration.  Their  presentation  surveys  discrete-time 
maintenance  models  and  continuous-time  maintenance  models.  In  discrete¬ 
time  maintenance  models,  a  component  is  monitored  and  a  decision  to  repair, 
replace  and/or  restock  the  component  is  made  at  discrete  intervals  in  time.  In 
continuous-time  maintenance  models,  actions  and  events  are  not  restricted,  a 
priori,  to  occur  only  within  a  discrete  subset  of  the  time  axis.  The  age 
replacement  preventive  maintenance  policy  falls  into  this  latter  category. 

Hence,  discussion  will  be  limited  to  this  area. 

Pierskalla  and  Voelker  [9]  note  that  research  has  extended  the  results  of 
the  Barlow  and  Proschan  [41]  optimal  age  replacement  model  for  the  infinite 
time  horizon  to  include  the  age  replacement  policy  under  a  total  discounted 
cost  criterion  [42]  as  well  as  including  an  age-dependent  cost  [43].  This  cost  is 
defined  as  the  increasing  burden  of  routine  maintenance  as  the  component 
ages  resulting  from  its  diminishing  productivity  or  reduced  salvage  value.  It  is 
also  noted  that  Glasser  [44]  has  extended  the  Barlow  and  Proschan's  optimal 
age  replacement  period,  T*,  for  the  truncated  normal,  gamma  and  Weibull 
density  functions.  Other  types  of  models  considered  include  age  dependent 
costs  that  are  incurred  at  discrete  times  [45],  and  a  generalization  to  include, 
c*,  a  maximum  cost  rate  that  when  exceeded  causes  the  component  to  be 
replaced  [46].  This  maximum  cost  rate,  c*,  is  set  by  the  decision  maker. 
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Multiple  component  systems  are  considered  by  Nakagawa  and  Osaki 
[47],  Their  research  considers  a  two-component  redundant  system  with 
identical  components.  One  component  operates  while  the  other  is  in  standby 
status.  The  component  in  standby  status  does  not  age.  Failure  of  the  system 
occurs  when  both  units  are  undergoing  maintenance  due  to  failure  or 
scheduled  preventive  action.  Under  the  conditions  where  mean  failure  repair 
time  is  greater  than  mean  preventive  maintenance  time,  they  find  the  optimal 
age,  T,  that  a  component  should  be  serviced.  Nakagawa  and  Osaki  find  the 
optimum  T  that  maximizes  the  limiting  availability.  Thus  an  infinite  time  horizon 
is  considered. 

Makabe  and  Morimura  [48,  49  and  50]  and  Morimura  [51]  consider  an 
age  replacement  model  where  the  component  is  replaced  at  the  failure. 

The  k-l  previous  failures  are  corrected  by  minimal  repair.  Under  minimal 
repair  a  component  is  repaired  to  operating  condition.  However,  the  repair  is 
assumed  to  not  effect  the  aging  of  the  component.  Hence,  after  the  minimal 
repair,  the  component  restarts  deterioration,  according  to  its  lifetime  density,  at 
the  same  point  in  time  that  it  began  minimal  repair.  For  components  with 
strictly  increasing  failure  rates,  the  optimal  policy  consists  of  t*  and  k.  If  the 
component’s  k^^  failure  occurs  before  t*  then  a  minimal  repair  is  completed  and 
replacement  is  completed  at  the  next  component  failure.  If  the  failure 
occurs  after  t*  then  a  replacement  is  carried  out.  Optimality  is  shown  for  the 
infinite  time  horizon  case  in  which  limiting  component  efficiency  (availability)  is 
maximized. 
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2.2.3  Preventive  Maintenance  Research  From  1976  to  Present 


In  "A  Survey  of  Preventive  Maintenance  Models  for  Stochastically 
Deteriorating  Single-Unit  Systems,"  Valdez-Flores  and  Feldman  [11]  (1989) 
present  a  state-of-the-art  survey  of  the  research  accomplished  in  preventive 
maintenance  mathematical  modeling  since  the  survey  accomplished  in  1975  by 
Pierskalla  and  Voelker  [9]  as  well  as  the  bibliographic  references  provided  by 
Osaki  and  Nakagawa  [52],  Sherif  and  Smith  [53]  and  Sherif  [54],  The 
maintenance  models  considered  in  this  survey  include  those  for  systems  that 
may  be  modeled  as  a  single  unit  or  component. 

Valdez-Flores  and  Feldman  classify  age  replacement  models  within 
what  they  term  "inspection  models".  This  group  of  models  considers  the  case 
where  the  status  of  the  component  or  system  may  not  be  continuously 
observed.  However,  it  may  be  possible  to  perform  an  inspection  to  ascertain 
the  status  of  the  system.  These  models  usually  assume  that  the  state  of  the 
system  is  unknown  unless  an  inspection  is  performed.  Initial  models  assumed 
that  the  inspection  information  was  perfect.  Later,  research  allowed  for  the 
case  where  inspections  were  not  perfect.  Thus,  age  replacement  preventive 
maintenance  research  has  progressed  now  to  include  inspection  policies. 
Optimal  policies  are  still  developed  by  minimizing  the  cost  rate  function  or  the 
expected  total  cost  per  inspection  cycle.  Barlow,  et  al  [37]  presented  this 
research  first. 
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Beichelt  [55]  extends  this  basic  model  by  finding  optimal  inspection 
times  for  the  cases  when  system  replacement  and  no  system  replacement  of  a 
failed  system  are  permitted.  He  also  obtains  optimal  inspection  schedules 
when  the  lifetime  distribution  of  the  system  is  partially  known  and  when  it  is 
completely  unknown  (system  mean  lifetime  is  considered  known).  Thus  this 
research  extends  the  inspection  model  to  cases  where  incomplete  or  unknown 
initial  conditions  prevail.  Luss  [56]  extends  this  further  to  the  case  where 
system  or  component  operational  state  is  multi-state  and  is  only  observable 
through  inspection.  The  model  is  generalized  for  L  states  of  component 
operational  condition.  Optimal  control  limit  policy  with  control  state  a  and  the 
optimal  inspection  intervals  for  states  0,  1, ...,  a  -  1  are  found  through  an 
iterative  procedure  by  assuming  exponentially  distributed  state  transition  times 
(Markovian  assumptions). 

Interestingly,  Rosenfield  [57  and  58]  examines  a  model  where 
component  deterioration  follows  a  discrete-time  Markov  chain.  The  operating 
and  replacement  costs  increase  with  state  number  and  inspections  are 
considered  perfect.  Kander  [59]  also  considers  inspection  for  a  component 
that  deteriorates  at  discrete  levels.  He  models  the  situation  using  semi-Markov 
processes  to  determine  the  optimal  inspection  schedule  that  minimizes  long- 
run  expected  cost  per  unit  time. 

Nakagawa  and  Yasui  [60]  present  an  algorithm  to  compute  sub-optimal 
inspection  times  for  components  that  do  not  have  exponential  failure  times. 
Their  research  demonstrates  the  validity  of  a  numerical  example  for  the  case 
where  the  component's  lifetime  distribution  is  Weibull.  A  recursive  approach 


25 


computes  successive  inspection  times  backwards  .  For  lifetime  hazard 
functions  characterized  by  an  increasing  hazard  rate,  Munford  [61  ]  shows  that 
an  inspection  policy  with  decreasing  intervals  between  successive  inspections 
as  a  function  of  component  age  are  better  than  strictly  periodic  inspection 
policies. 

Minimal  repair  models  have  also  been  considered.  Park  [62]  presents 
an  entirely  new  concept  of  age  replacement  under  minimal  repair  as  that 
introduced  by  Barlow  and  Proschan  [13].  Park  [62]  finds  the  number  of  failures 
and  minimal  repairs  before  the  system  is  replaced  instead  of  a  fixed  age  to 
replacement.  He  develops  a  closed-form  solution  of  his  model  for  the  case 
where  the  system's  lifetime  distribution  is  Weibull.  Numerically,  Park  [62] 
shows  his  policy  provides  a  lower  long-run  expected  cost  per  unit  time  than  a 
fixed-age  replacement  policy.  However,  a  mathematical  proof  is  not  offered. 

Recent  research  includes  Nachlas,  1989  [63]  and  Lam  and  Yeh,  1994 
[64].  Nachlas  develops  a  constrained  optimization  model  for  selecting 
replacement  policies  for  both  age  and  block  replacement.  This  model's 
solution  finds  the  replacement  policy  that  minimizes  the  cost  rate  function 
constrained  by  a  confidence  interval  on  required  availability.  Nachlas's  model 
assumes  failure  repair  and  preventive  maintenance  times  are  random  but 
identically  and  independently  distributed.  Lam  and  Yeh  present  algorithms  for 
deriving  optimal  maintenance  policies  to  minimize  the  mean  long-run  cost-rate 
for  continuous-time  Markov  deteriorating  systems.  The  model  assumes  the 
deterioration  state  is  only  known  through  inspection.  They  consider  five 
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maintenance  strategies.  These  are  failure  replacement,  age  replacement, 
sequential  inspection,  periodic  inspection  and  continuous  inspection. 


2.2.4  Summary  of  Preventive  Maintenance  Research 

The  models  described  above  consider  the  infinite  time  horizon  case  in 
which  cost  rate  functions  are  developed  and  minimized.  The  alternative 
approach  presented  uses  Markov  chain  analysis.  In  the  case  of  the  infinite 
time  horizon  case,  cost  and/or  availability  functions  with  respect  to  time  are  not 
developed.  Under  Markovian  assumptions,  the  lifetime,  failure  repair  and 
preventive  maintenance  probability  density  functions  are  exponential.  Under 
these  restrictive  assumptions,  cost  and  availability  functions  with  respect  to 
time  may  be  developed.  Applicable  density  functions  will  be  discussed  later  in 
this  review.  However,  it  is  evident  that  the  instantaneous  availability  function, 
A{t),  is  a  prime  research  area  for  the  cases  where  lifetime,  repair  and 
preventive  maintenance  distributions  are  not  all  exponentially  distributed. 


2.3  Analytical  Techniques 

The  summary  of  the  preventive  maintenance  models  presented  above 
shows  that  with  the  exception  of  Markov  models,  the  optimal  preventive 
maintenance  policy  is  derived  by  developing  an  expected  total  cost  rate 
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function.  This  function  is  minimized  to  find  the  optimal  policy.  The 
development  of  this  function  requires  that  an  infinite  time  horizon  be  assumed. 
Often,  this  is  not  realistic  as  system  lifecycles  consist  of  a  finite  operational 
period.  This  section  summarizes  analytical  techniques  and  methods  pertinent 
to  developing  the  component  instantaneous  availability,  A(t),  for  an  age 
replacement  policy. 


2.3.1  Renewal  Theory 

Renewal  theory  deals  with  stochastic  processes  that  have  renewal 
points.  At  these  renewal  points,  the  stochastic  process  probabilistically  starts 
again  or  renews  itself.  A  renewal  process  is  therefore  made  up  of  a  finite 
number  of  random  periods,  n,  each  with  density  function,  /^(f)  Vi  =  1, n, 
which  reoccur  in  a  specific  order  (cycle)  as  the  process  progresses.  Thus  the 
renewal  process  consists  of  these  well  defined  cycles  that  are  probabilistically 
identical.  Interestingly,  renewal  theory  began  as  the  study  of  probability 
problems  connected  with  the  failure  and  replacement  of  components,  such  as 
light  bulbs  [26].  Barlow  and  Proschan  [1 3]  and  Cox  [26]  offer  a  mathematical 
development  of  renewal  theory. 

One  important  property  stems  from  the  definition  of  a  renewal  process 
given  above.  Since  a  time-based  renewal  process  deals  with  some  finite 
number  of  random  time  periods  that  occur  in  a  specific  order  during  a  cycle,  it 
is  often  necessary  to  find  the  density  function  for  the  time  between  renewals. 
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This  may  be  accomplished  by  taking  the  convolution  of  the  probability  density 
functions  represented  in  a  renewal  cycle,  /,(f)  Vi  =  1, n.  For  the  time-based 
renewal  process,  this  is  known  as  the  density  function  for  time  between 
renewals  and  is  designated  h{t).  Mathematically,  h{t)  is  stated  as 

h{t)  =  /,(f)  (8)  /2(f)  (8) ...  (8)  <8>  Ut)  .  (2.1 ) 

where  the  symbol,  (8 ,  represents  the  convolution  operator.  The  convolution  of 
two  functions,  /j(f)  ®  fo{t),  is  mathematically  defined  as 

fi(t)  (8)  /2(f)  =  Jofiit)  f2{t  -u)du  =  /0/2  W  flit  -u)du. 


2.3. 1.1  Laplace  Transforms 

The  probability  density  function,  h{t),  may  be  found  mathematically  by 
taking  the  convolution  as  defined  in  equation  2.1. 

The  Laplace  transform  convolution  theorem  provides  a  method  for 
finding  h{t)  without  explicitly  evaluating  the  multiple  integrals  required  for 
equation  2.1 .  If  we  designate  the  Laplace  transform  of  a  function,  £[/(f)],  as 
/*(s)then  the  following  statements  hold  [97]: 

(1)  ^[fit)\  =  r{s) 

(2)  £-'|r(3)|=/(i) 
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(3)  Mt)^m  =  z-^[f;{s)f;is)]. 


For  many  cases,  the  inversion  of  the  Laplace  transform  into  a  time 
domain  function  may  be  found  through  tables,  partial  fraction  decomposition, 
contour  integral  evaluation  as  well  as  other  methods.  However,  the  case  may 
arise  where  the  inverse  of  h*{s)  may  not  be  found  exactly.  In  this  case, 
approximation  methods  must  be  used.  These  methods  include  numerical 
inversion  techniques  evaluated  by  Davies  and  Martin  [65]  as  well  as  those 
evaluated  and  implemented  in  Mathematica  (Wolfram  Research  Inc.)  by 
Cheng,  et  al  [66].  Several  numerical  inversion  methods  for  Laplace  transforms 
are  evaluated.  These  include  the  Stehfest  method  [67],  the  Papoulis  legendre 
polynomial  method  [68],  the  Durbin  method  [69],  the  Crump  method  [70],  the 
Weeks  method  [7 1  ],  and  the  Piessens  gaussian  quadrature  method  [72  and 
73].  These  methods  provide  an  estimate  of  the  time  domain  function  at  a 
specific  value  of  time,  f{t  -  c).  Functions  using  the  Laplace  transform 
function,  /*(s),  are  developed  and  implemented  [66]. 


2.3.2  Lifetime  Probability  Density  Functions 

A  component's  lifetime  probability  density  function,  defines  the 
probability  density  function  for  time  to  component  failure.  This  function  defines 
the  probability  of  component  failure  at  time  t  +  e,  where  e  is  an  infinitesimally 
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small  increment  in  time.  Much  research  has  been  accomplished  concerning 
the  proper  probabilistic  modeling  of  the  component  lifetime  density  function. 
The  most  frequently  used  density  function  is  the  exponential  density  function. 
This  density  function  seems  to  model  electronics  well.  The  most 
advantageous  aspect  of  this  probability  density  function  is  its  mathematical 
tractability  and  it  is  implemented  often  for  this  reason  alone  [88].  Justification 
for  the  use  of  the  exponential  lifetime  density  may  be  found  in  [74,  75,  76,  77, 
78,  79,  80,  81  and  82]. 

Other  references  such  as  Bell,  et  al  [83]  maintain  that  density  functions 
such  as  the  truncated  normal,  log-normal  and  Weibull  density  functions  are 
applicable  in  many  cases.  Bell,  et  al  [83]  claim  that  a  significant  portion  of 
aircraft  parts  follow  a  normal  lifetime  density  function.  Lie,  et  al  [10]  make  the 
following  statements  concerning  lifetime  density  functions; 

(1 )  In  many  practical  situations,  lifetime  density  functions  are  markedly 
skewed  and  thus  are  not  normally  distributed. 

(2)  The  gamma  and  Weibull  density  functions  provide  the  skewed 
characteristic  discussed  above. 

(3)  The  Weibull  density  function  is  known  to  be  suitable  in  describing 
fatigue  failure  such  as  that  occurring  in  vacuum  tubes  and  ball-bearings. 

(4)  The  log-normal  density  function  seems  to  fit  repair  times  than  failure 

times. 

The  Weibull  density  function  is  justified  for  use  as  a  lifetime  density 
function  in  [83,  84,  85  and  86]  and  the  truncated  normal  in  [83].  Lie,  et  al  [10] 


31 


provide  additional  references  on  other  lifetime  density  functions  such  as  the 
Erlang,  gamma,  Rayleigh,  uniform,  extreme  value  and  arbitrary  functions. 

Proper  modeling  of  a  component's  lifetime  density  must  consider  the 
component's  family  (mechanical,  electronic,  etc.)  characteristics  and  the  failure 
mechanisms  at  work  on  and/or  within  the  component. 


2.3.3  Repair  Time  Probability  Density  Functions 

Repair  time  density  functions  have  usually  been  approximated  by  an 
exponential  density  function  for  analytical  purposes  [87].  Rohn  [88]  considers 
the  repair  times  of  complex  electronic  equipment.  He  maintains,  in  the  case  of 
electronic  repair,  that  the  repair  times  consist  of  a  high  frequency  of  short  times 
and  a  few  long  repair  times.  This  characteristic,  according  to  Rohn,  suggests 
an  exponential  lifetime  density  function. 

Other  studies  concerning  airborne  radar  equipment  and  ground 
equipment  for  surface-to-air  missile  systems  have  shown  that  the  log-normal 
density  function  fits  the  observed  repair  times  best  [89  and  90].  Additional 
references  justifying  the  use  of  the  log-normal  probability  density  function  to 
model  repair  times  include  [91,  92,  93  and  94], 
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CHAPTER  3 


PROBLEM  STATEMENT 


3.1  Stochastic  Process  Description 


Failure  Cycle  Fk'eventive  Maintenance  Cycle 


Operational  Period 

Repair 

Operational  Period 

¥M. 

••• 

t=T 

FailurebefareHiiieT  Unit  survives  iindl Tine  T 


Figure  3.1  Age  Replacement  Preventive  Maintenance  Stochastic  Process 

The  problem  analyzed  in  this  research  is  depicted  in  Figure  3. 1 .  This 
figure  represents  a  stochastic  process  based  upon  an  age  replacement 
preventive  maintenance  policy  for  a  continuously  demanded  component.  In 
this  situation,  the  component  is  repaired  or  replaced  if  it  fails  before  T  time 
units  have  passed  since  the  last  maintenance  action.  This  event  is  called  a 
failure  cycle.  The  component  is  preventively  maintained  whenever  T  time  units 
have  passed  since  the  last  maintenance  action.  This  event  is  called  a 
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preventive  maintenance  cycle.  The  stochastic  process  consists  of  randomly 
occurring  failure  cycles  and  preventive  maintenance  cycles.  Figure  3.1  does 
not  imply  an  ordering  in  the  occurrence  of  failure  and  preventive  maintenance 
cycles.  The  occurrence  of  these  cycles  is  completely  random. 

Failure  repairs  and  preventive  maintenance  actions  are  perfect;  thus  at 
the  completion  of  a  failure  repair  or  preventive  maintenance  action,  the 
component  is  considered  to  be  as  good  as  new.  At  these  points  in  time  we 
have  component  renewal  and  thus  have  a  stochastic  process  exhibiting 
renewal  points.  At  these  points  the  stochastic  process  renews  itself  and 
probabilistically  repeats. 

The  notation  used  in  this  research  is  shown  below.  The  symbol,  =  ,  is 
read  as  "is  defined  as". 

Notation: 

/jr  (f)  =  known  component  lifetime  density 

=  known  component  lifetime  distribution,  /o/iCf)  di 

RL{t)  =  known  component  survivor  function,  1  - 

g^{t)  =  known  repair  time  density  after  component  failure 
=  known  preventive  maintenance  time  density 

T  =  age  replacement  preventive  maintenance  period. 

The  stochastic  processes  governing  the  component's  lifetime,  the  time 
required  to  effect  repair  after  a  component  failure  and  the  time  required  to 
carry  out  a  preventive  maintenance  action  are  considered  to  be  governed  by 
known  probability  density  functions,  which  are  all  unique  and  independent 
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functions.  Furthermore,  it  is  assumed  these  probability  density  functions  are 
integrable  in  closed  form  to  obtain  their  respective  distribution  functions. 


3.2  Problem  Discussion 


As  stated  earlier,  past  research  has  concentrated  on  finding  an  optimal 

replacement  period,  T^,  for  component  use  over  an  infinite  time  horizon.  This 

optimal  replacement  period  maximizes  the  limiting  availability  of  the 

component.  In  this  research,  the  availability  function,  A{t),  is  defined  to  be  the 

probability  the  component  operates  at  time  t.  The  limiting  availability,  A,  is 
defined  to  be  lim  A{t). 

t—*oo 

An  infinite  time  horizon  for  the  component's  expected  life  cycle  is  not 
realistic.  Realistically,  a  component  is  used  over  a  prescribed  time  period.  For 
a  finite  component  life  cycle,  the  average  availability  over  the  time  period  may 
be  a  more  appropriate  measure  to  be  maximized.  The  average  availability, 
Aavgir),  is  defined  as 


A^,{T)  =  {l)S;A{t)dt.  (3.1) 

Finding  a  value  for  the  age  replacement  period,  T,  that  maximizes  Aavgir) 
requires  an  expression  for  the  availability  function,  A{t). 
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3.3  Development  of  the  Availability  Function 

A  notable  feature  of  the  stochastic  process  under  consideration  is  the 
unique  probability  density  functions  representing  the  time  required  for  repair 
(after  component  failure)  and  the  time  required  for  preventive  maintenance 
after  T  time  units  have  passed  since  the  last  maintenance  action.  Previous 
results  only  consider  the  expected  values  of  these  times  since  the  analysis  is 
for  an  infinite  time  horizon  (t-*oo).  For  finite  time  horizons,  the  nature  of  the 

probability  density  functions  governing  these  respective  service  times  must  be 
taken  into  account. 

To  accommodate  this,  a  renewal  theory  approach  is  taken  to  find  the 
availability  function  for  the  specified  stochastic  process.  This  approach  is 
based  upon  Laplace  transform  theory  as  well  as  the  application  of  conditional 
probability  theory.  The  result  of  the  initial  analysis  is  the  Laplace  transform  of 
the  availability  function,  A*{s),  which  is  defined  as 

A*(s)  =  J^(e-^^)A(t)dt.  (3.2) 

The  availability  function,  y4(f),  may  be  found  by  taking  the  inverse  Laplace 
transform  of  A*(s).  It  is  important  to  realize  that  throughout  this  development, 
T,  the  age  replacement  preventive  maintenance  period,  is  treated  as  a 
parameter.  After  inverting  A*(s),  we  obtain  the  availability  as  a  function  of  both 
time,  t,  and  the  age  replacement  preventive  maintenance  period,  T. 
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As  shown  in  Figure  3.1,  two  types  of  sample  paths  occur.  The  first  type 
of  cycle  shown  is  the  "failure  cycle".  In  this  event  the  component  fails  before  T 
time  units  have  elapsed  since  the  last  repair  or  preventive  maintenance  action. 
At  failure,  a  failure  service  action  immediately  begins.  The  component's 
operational  or  uptime  is  always  less  than  T  in  this  cycle.  Thus,  the 
component's  lifetime  probability  density  function  is  truncated  at  time  T. 

The  second  type  of  cycle  shown  is  the  "preventive  maintenance  cycle". 
In  this  event  the  component  survives  exactly  T  time  units  since  the  last  repair 
or  preventive  maintenance  action.  After  the  component  reaches  the  "age"  of  T 
the  component  undergoes  a  random  repair  period  (preventive  maintenance 
action).  Thus  the  component's  operational  period  or  uptime  in  this  case  is  a 
constant  T  time  units. 

The  service  time  periods  are  not  truncated  and  may  range  on  [0,oo].  In 
the  failure  cycle,  we  are  dealing  with  an  operating  period  of  f  <  T  and  in  the 
preventive  maintenance  cycle  we  have  a  constant  operating  period  of  f  =  T. 
These  limitations  introduce  special  cases  of  the  Laplace  transform  that 
complicate  our  expression  for  the  Laplace  transform  of  the  availability  function, 
A’(s).  The  following  sections  review  pertinent  properties  of  the  Laplace 
transform  applicable  to  analyzing  the  age  replacement  preventive  maintenance 
policy. 
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3.4  Laplace  Transforms 


Equation  (3.2)  gives  the  classical  definition  of  the  Laplace  transform.  In 
some  cases  we  may  experience  a  truncated  time  period.  Specifically,  the  time 
domain  variable  may  have  a  limited  range.  This  is  true  for  the  stochastic 
process  under  investigation  and  discussed  above.  For  example  consider  the 
"failure  cycle"  and  the  component's  lifetime  probability  density.  Let  f^(t)  be  the 
non-truncated  probability  density  function  for  the  component's  lifetime  in  which 
t  has  range  [0,oo].  If  operating  period,  t,  is  now  limited  to  the  range,  [0,T],  then 
the  truncated  lifetime  probability  density  function  is 

0  otherwise  . 

The  partial  or  truncated  Laplace  transform  for  f{t)  is 

r(«.r)  =  . 

Basic  Laplace  transform  properties  include  the  useful  convolution 
theorem.  The  convolution  of  two  functions,  f{t)  and  g{t)  is 

fit)  ®  git)  =  Jo  fit)  git  -u)du  =  Jo  f{t  -  u)  g{t )  du , 
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where  the  symbol,  ® ,  represents  the  cxjnvolution  operator.  The  Laplace 
transform  convolution  theorem  provides  the  following  relationship: 

fit)  ®  git)  =  p’(s)] . 

where  £'''[•]  is  the  inverse  Laplace  transform  of  the  argument  within  the 
brackets  [97], 

Further  discussion  of  the  partial  Laplace  transform  is  appropriate.  Due 
to  the  limited  operating  time  range  of  the  component,  notably ,  0  <  f  <  T.  we 
are  forced  to  reconsider  the  limits  on  the  Laplace  transform  for  the  truncated 
lifetime  density.  This  range  limitation  produces  the  result  shown  for  /*(s,  T) 
above.  Section  4.3  produces  a  significant  result  for  the  partial  Laplace 
transform. 


3.4.1  ADProximate  Laplace  Transforms 

Many  of  the  most  useful  lifetime  probability  distributions  as  well  as 
service  time  distributions  do  not  have  exact  Laplace  transforms.  An  approach 
similar  to  the  one  used  to  approximate  Laplace  transform  inversion  is  often 
used.  The  time  domain  function,  f{t),  is  expanded  into  a  series  of  functions 
with  known,  exact  Laplace  transforms.  Thus  an  approximate  Laplace 
transform,  /*(s),  may  be  found. 
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If  the  series  of  the  approximate  Laplace  transform  contains  many  terms, 
then  the  complexity  of  the  A*{s)  increases.  Thus,  there  is  a  trade-off  between 
the  producing  an  accurate  estimated  Laplace  transform  and  obtaining  a 
manageable  expression  for  A*  (5).  Issues  pertaining  to  the  inversion  of  A*(s)to 
obtain  the  time  domain  function,  A{t)  are  discussed  next. 


3.4.2  Approximate  Laplace  Transform  Inversion 

Due  to  the  complicated  nature  of  the  Laplace  transform  of  the 
availability  function,  4*(5),  an  exact  inversion  may  not  exist  or  may  not  be 
easily  obtainable.  In  these  cases,  it  is  necessary  to  construct  approximations 
of  to  find  an  approximate  expression  for  A{t).  Numerical  approximation 
methods  have  been  developed.  These  methods  are  applied  through  the 
assistance  of  Mathematica  [96]  to  obtain  the  approximate  Laplace  transform 
inversion  of 

Conceptually  these  methods  attempt  to  expand  the  non-invertible  /‘(s) 
into  functions  of  s  that  have  exact  transforms.  Various  methods  exist  for  this 
expansion  and  these  methods  are  evaluated  for  use  in  this  analysis. 
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3.5  Maximizing  Average  Availability  Over  A  Finite  Time  Horizon 


After  obtaining  an  approximate  expression  for  A{t),  the  next  step 
attempts  to  obtain  a  value  for  T  that  maximizes  the  Aavg{t)  through  equation 
(3.1).  Two  options  are  investigated.  The  first  step  involves  inverting,  i4*(s). 
the  Laplace  transform  of  the  availability  function.  Recall  that  equation  (3.1 ) 
defines  Aavg{t)  and  involves  /q  A{u)  du.  A  useful  property  of  the  Laplace 
transform  allows  the  direct  development  of  Aavg{t)  through  a  simple 
transformation  of  >l*(s)  before  inversion.  If  f{t)  is  a  function  with 
fg  fit)  dt  —  F(t)  and  Laplace  transform,  /’(s),  then  the  following  holds: 

£■'[  ;  rWI  =  F(t) . 

Once  an  expression  is  obtained  for  .A* (5),  we  simply  multiply  by  (i)  and  invert 

to  obtain  A{t)  dt  [97].  Finally,  to  obtain  Aavg{t)  we  multiply  the  resulting 
function,  representing  A(t)  dt  by  (i).  This  result  gives  a  functional 

representation  of  Aavg{t).  Once  this  function  is  obtained,  optimization 
operations  may  be  performed  to  find  values  of  T  that  maximize  Aavg{t)  for 
specific  values  of  t. 

However,  the  resulting  function  obtained  for  Aavg{t)  may  not  be 
differentiable  or  of  a  form  that  is  convenient  to  manipulate.  An  alternative 
method  is  to  perform  a  numerical  search.  In  this  case,  we  would  select  a 
specific  time  horizon,  t,  and  perform  the  integration  for  various  values  of  T.  (It 
is  assumed  specific  values  have  been  given  to  all  probability  density  function 
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parameters.)  Ranges  for  the  possible  values  of  T  are  first  limited  by  the  time 
horizon,  t.  It  is  evident  that  1  <t,  because,  if  T  >  t,  then  no  preventive 
maintenance  would  be  performed  during  the  time  horizon  period.  Also  the 
value  for  T  obtained  for  the  infinite  time  horizon  case  should  provide  some  limit 
on  T  for  the  finite  time  horizon  case.  However,  it  is  not  clear  whether  this  may 
be  an  upper  or  lower  bound. 

Clearly,  the  most  interesting  cases  include  the  condition  where  aging  is 
exhibited  by  the  component.  In  reliability  theory,  this  is  called  increasing 

hazard.  The  hazard  function,  z{t)  =  is  strictly  increasing  if  >  0.  If 

the  component  exhibits  constant  hazard  or  decreasing  hazard,  preventive 

maintenance  will  adversely  effect  the  component's  availability.  Under  constant 
hazard  (^^=0),  no  matter  how  long  the  component  has  survived,  it  is  always 

as  good  as  new  (exponential  lifetime  density  function).  Thus  there  would 
never  be  an  age  at  which  the  component  is  purposely  brought  out  of  operation 
for  a  preventive  maintenance  action.  Under  strictly  decreasing  hazard 

<  0),  the  longer  the  component  survives  the  more  reliable  the  component 

becomes.  Thus,  it  is  easily  seen  that  if  the  component  actually  becomes  more 
reliable  with  time,  then  preventive  maintenance  would  adversely  affect 
component  availability.  If  we  consider  how  the  probability  of  failure  changes 
with  time,  we  will  find  it  actually  decreases  the  longer  the  component  survives. 
Some  references  term  this  probability  of  failure  the  "force  of  mortality"  or  "age- 
specific  failure  rate"  [26].  For  components  with  strictly  increasing  hazard,  the 
force  of  mortality  increases  as  the  component  ages;  for  components  with 
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constant  hazard  the  force  of  mortality  remains  constant  as  the  component 
ages;  for  components  with  strictly  decreasing  hazard,  the  force  of  mortality 
decreases  as  the  component  ages. 

Additionally,  the  stochastic  nature  of  the  failure  repair  time  and  the 
preventive  maintenance  action  time  must  be  considered.  In  this  analysis,  we 
consider  the  case  in  which  the  failure  repair  time  distribution  is  stochastically 
greater  than  the  preventive  maintenance  time  distribution.  A  random  variable 
X  is  stochastically  greater  than  another  random  variable  Y  if  the  following  is 
true: 

Prob[X  >  z]  >  Prob\Y  >  z],  f  or  all  z, 

or  equivalently, 

Fx{z)  <  Gy{z),  for  all  z, 

where  Fxiu)  and  Gy  (u)  are  the  cumulative  density  functions  for  the  random 
variables,  X  and  Y,  respectively  (see  p.78  of  [95]).  Logically,  it  is  easily 
argued  that  if  this  condition  does  not  exist,  then  preventive  maintenance  will 
adversely  affect  the  availability  of  the  component.  If  the  preventive 
maintenance  action  time  distribution  were  stochastically  greater  than  or  equal 
to  the  failure  repair  time  distribution,  the  option  would  be  to  let  T-»oo,  and  thus 
perform  no  preventive  maintenance.  In  this  case,  the  component  is  only 
repaired/replaced  upon  failure.  Thus  the  important  conditions  to  be  analyzed 
include  component  lifetime  distributions  with  strictly  increasing  hazard  and 
repair  time  distributions  that  are  stochastically  greater  than  preventive 
maintenance  action  time  distributions. 
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Another  point  needs  to  be  addressed  concerning  repair  and  preventive 
maintenance  time  distributions.  Only  distributions  from  the  same  families  are 
considered  here.  For  example,  if  the  failure  repair  distribution  is  modeled  as 
exponential,  then  the  preventive  maintenance  time  distribution  is  also  modeled 
as  exponential.  Analysis  considering  distributions  with  shape  parameters  (i.e., 
Weibull)  use  the  same  shape  parameter  for  both  distributions.  In  these  cases, 
only  distributions  with  strictly  increasing  hazard  are  considered.  Repair  time 
distributions  exhibiting  these  characteristics  are  interpreted  as  follows:  given 
that  a  repair  task  has  taken  t  time  units,  the  probability  of  finishing  the  task  in 
the  next  instant  of  time  increases  as  t  increases. 


3.6  Summary  of  Analytical  Approach 

The  stochastic  process  shown  in  Figure  3.1  is  analyzed  as  follows. 

First,  a  renewal  theory  based  model  is  developed.  The  objective  of  this  model 
is  to  evaluate  the  availability  of  the  system  or  component.  Expressions  for  the 
availability  of  the  component  or  system  are  developed.  These  expressions 
require  the  evaluation  of  the  process  renewal  density,  and/or  the 
probability  density  function  for  the  time  between  renewals,  h{t).  The  Laplace 
transform  is  used  to  find  the  convolutions  of  the  probability  density  functions. 

Once  the  Laplace  transform  of  the  availability  function,  A*(s),  has  been 
developed,  it  must  be  inverted.  Numerical  inversion  packages  are  investigated 
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for  application  to  this  analysis.  Also  methods  to  invert  functions  not  having  an 
exact  closed-form  Laplace  transform  are  investigated. 

These  Laplace  transform  and  inversion  methodologies  are  applied  to 
specific  cases.  These  cases  include  an  all  exponential  case  for  validation  and 
a  Weibull  lifetime  density  with  exponentially  distributed  repair  times  and 
preventive  maintenance  times.  In  all  cases,  the  distributions  are  unique  and 
independent. 

The  following  ratios  are  established  for  this  analysis: 

_  _ Mean  Lifetime _ 

^  Mean  Failure  Repair  Time 

^  M ean  Failure  Repair  T ime 

MeanP.M.  ActionTime 

Various  values  of  p  and  6  are  placed  in  the  model  and  evaluated.  Typical 
values  of  p  range  from  two  up  to  ten  thousand  while  values  of  6  range  from  two 
up  to  ten.  Lower  values  of  p  are  expected  to  provide  more  noticeable  changes 
in  the  average  availability  as  T  is  changed.  This  is  easily  understood  by  noting 
that  the  lower  the  value  of  p,  the  longer  the  repair  and  preventive  maintenance 
periods  are  in  relation  to  the  operative  lifetime.  Thus  we  would  expect  the 
average  availability  to  be  less  for  lower  values  of  p  given  all  other  parameters 
are  held  constant.  The  values  obtained  for  T  for  finite  time  horizons  are 
compared  to  the  values  obtained  for  T  for  infinite  time  horizons. 
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The  infinite  time  horizon  values  also  allow  the  model  to  be  evaluated 
and  partially  validated  for  large  values  of  t.  Since  the  lim  Aavg{t)  =  Urn  A{t), 

t—*oo  i— >oo 

we  may  take  a  large  time  horizon  and  evaluate  the  average  availability  at  the 
value  of  T  found  from  the  infinite  time  horizon  and  compare  it  to  the  limiting 
availability  obtained  from  the  infinite  time  horizon  model.  The  model  is  also 
checked  at  very  small  values  of  time  where  the  availability  of  the  model  should 
be  close  to  unity. 
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CHAPTER  4 


MODEL  DEVELOPMENT 


4. 1  Component  Availability 

Barlow  and  Proschan  [7]  define  the  availability  function,  A{t),  as  the 
probability  that  a  component  is  functioning  at  time  t.  They  develop  the  integral 
equation  for  the  availability  function  of  an  entity  modeled  by  an  alternating 
renewal  process.  This  process  is  the  repair  upon  component  failure  case  in 
which  no  preventive  maintenance  is  performed  (T-»oo).  The  availability 
equation  for  the  repair  upon  failure  (alternating  renewal)  process  is 

A(t)  =  Riit)  +  Jq  Ri(t  -  u)  m(u)  du  ,  (4.1 ) 

where  Riit)  is  the  component's  survivor  function  (1-Fi(f))  and  m{t)  is  the 
renewal  density  function  [7].  This  expression  for  the  availability  as  a  function  of 
time,  t,  may  be  interpreted  as  follows;  the  probability  the  component  is 
functioning  at  time  t  is  equal  to  the  probability  that  the  unit  has  not  failed  since 
the  time  the  component  was  put  into  service  (represented  by  the  Ri{t)  term) 
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plus  the  probability  that  a  renewal  occurred  at  time  u,  (represented  by  the  m(u) 
term)  and  no  component  failure  occurred  in  the  remaining  time  period 
(represented  by  the  Riit  -  u)  term).  The  integral  sums  the  probabilities  of 
these  possible  events  over  the  time  interval  [0,t]. 


4.1.1  Laplace  T ransform  of  the  Availability  Integral  Equation  for  the 
Repair  Upon  Failure  Case 

Equation  (4.1)  gives  an  integral  expression  for  finding  the  availability  of 
a  component  that  is  repaired  upon  failure.  This  equation  may  be  solved  by 
taking  the  Laplace  transform  of  both  sides.  Before  proceeding,  recognize  that 
the  second  term  on  the  right-hand  side  of  equation  (4.1)  is  a  convolution  of  the 
component  survivor  function  and  the  renewal  density  function  of  the  process. 
This  convolution  may  be  found  by  taking  the  inverse  Laplace  transform  of  the 
product  of  the  Laplace  transforms  of  the  survivor  function  and  the  renewal 
density  (Tierney,  p.251,[97]).  This  relationship  is 

(8)m(f)  =  Jq  R{t  —  u)  m{u)  du  —  £'^[Ri{s)  m*{s)], 

where /*(s)  =  £[f{t)]  =  Using  this  relationship  and  taking  the 

Laplace  transform  of  both  sides  of  equation  (4.1)  we  obtain 

A*(s)  =  Rl{e)  +  {Rl{s)m*{3))  . 
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Factoring  the  RKs),  we  obtain 


=  Rl{s){l+m*{s)).  (4.2) 

The  renewal  density,  m(t)  is  defined  by  Cox  [26]  and  Barlow  and  Proschan  [7] 
as  the  sum  of  the  n-fold  convolutions  of  the  probability  density  function  for  the 
time  between  renewals,  h{t).  Proofs  may  be  found  in  the  references  above  as 
well  as  in  Ross,  p.  297  [101],  This  relationship  is 

oo 

m{t)  =  ,  (4.3) 

n=l 


where  {/i(f)}"  is  the  n-fold  convolution  of  h{t).  By  taking  the  Laplace  transform 
of  both  sides  of  equation  (4.3) ,  Cox  [26](p.54)  has  shown  a  simplification  of 
this  infinite  sum  as 


m-(3)  =  Z{h'{:)r  . 

n=l 


(4.4) 


where  in  this  case  (/i*(s))"  is  h*{s)  raised  to  the  n—  power.  Barlow  and 
Proschan  [7]  (p.166)  reduce  equation  (4.4)  to 


Tn*{s) 


h'js) 

l-ft’(s) 


(4.5) 
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Substituting  this  result  into  equation  (4.2),  the  following  is  obtained. 

Simplifying,  the  desired  result  is 


Recall  that  i?2(s)  is  the  Laplace  transform  of  the  component's  survivor  function 
and  that  h*(s)  is  the  Laplace  transform  of  the  probability  density  function  for  the 
time  between  renewals  in  the  process.  The  main  significance  and  advantage 
of  equation  (4.6)  is  that  the  Laplace  transform  expression  for  the  component 
availability  function  does  not  require  the  Laplace  transform  of  the  renewal 
density  function,  m*(s). 


4.2  The  Availability  Function  for  a  Component  Age  ReplanemAnt 
Preventive  Maintenance  Pniiry 

This  development  begins  by  breaking  the  time  scale  into  two  intervals. 
The  first  includes  the  interval  where  the  availability  function  of  the  component 
is  desired  for  t  <  T,  where  t  is  the  time  variable  and  T  is  the  age  replacement 
period.  Recall  that  the  component  is  preventively  maintained  whenever  T  time 
units  pass  since  the  last  maintenance  action,  whether  failure  repair  or 
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preventive  maintenance  was  last  performed  The  second  interval  considers 
t>T. 

For  f  <  T,  a  component  repair/replacement  due  to  preventive 
maintenance  will  not  have  occurred.  All  repair/replacements  during  this  time 
interval  are  due  to  component  failure.  Thus,  for  <  <  T,  a  repair  upon  failure 
process  (alternating  renewal)  is  experienced  and  the  availability  function  as  in 
equation  (4. 1 )  is  expressed  as 

A(t)  =  Riit)  +  /o  Riit  —  w)  m{u)  du ,  for  f  <  T. 

fort  >  T,  component  renewal  may  occur  due  to  component  failure  or  by 
preventive  maintenance.  In  this  case,  the  farthest  point  back  in  time  a  renewal 
may  occur  from  time  f,  is  f  -  T  time  units.  This  condition  holds  since  an  age 
replacement  policy  of  T  time  units  is  in  place.  Thus,  our  time  interval  is  limited 
to  [t  -  T,  f].  The  expression  for  the  availability  function  for  f  >  T  is 

■^(t)  =  ft-j  —  u)  m{u)  du,  for  f  >  T.  (4.7) 

Equation  (4.7)  may  be  interpreted  as  follows:  for  f  >  T,  the  probability  that  the 
component  functions  at  time  t,  A{t),  is  equal  to  the  probability  that  a  renewal 
occurs  at  time  u  (represented  by  Tn(tx))  and  the  component  survives  the 
remaining  period  of  time  (represented  by  Ri{t-u)).  The  renewals,  which 
occur  due  to  component  failure  or  scheduled  preventive  maintenance,  may 
occur  continuously  in  time  during  the  interval  [t  -T,t].  Only  this  time  interval 
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is  of  interest  in  the  case  f  >  T,  since  the  longest  time  that  may  elapse  after  a 
renewal  is  T  time  units.  This  particular  event  has  a  renewal  (preventive 
maintenance  in  this  case)  occurring  at  time  t-J  and  no  failures  in  the 
remaining  period  of  time. 

Summarizing,  the  availability  function  for  a  component  age  replacement 
preventive  maintenance  policy  is 


A{t)  = 


fo  m(u)  du,  for  t  <  T, 


j  RL(t-u)Tn{u)du,  fort>T. 


(4.8) 


4.3  Laplace  Transform  for  the  Availability  Function  for  a  Component 
Age  Replacement  Preventive  Maintenance  Policy 

The  following  development  shows  the  construction  of  the  Laplace 
transform  of  the  availability  function,  A(t),  stated  in  equation  (4.8).  The  proof 
is  started  by  noting  the  definition  of  the  Laplace  transform  for  the  availability 
function  is 


A*{s)  =  A{t)dt. 
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Applying  (4.8)  to  this  relationship,  the  following  is  obtained: 


=  fj  e""'  Ri(t)  dt  +  fj (/o.Ri(f  -  u)  m(u)  du^  dt 

+  ‘  ^ Ri{t  —  u)  m{u)  du^ dt 

=  Jq  Ri(t)  dt  +  fJ Jq  Ri(t  —  u)m(u)du  dt 

+  J^_j  Riii  ~  '^)  Tn{u)  du  dt  . 

Changing  the  order  of  integration  on  the  last  two  terms  involving  the  double 
integrals,  the  following  is  obtained: 


A*(s)  =  Qe  Ri{t)dt  +  fjfje  Ri(t  — u)m(u)dt  du 

+  Ri(t  —  u)  m{u)  dtdu  +  Ri{t  —  u)  m{u)  dtdu 

=  fJ Ri(t)  dt  +  fjm(u)  Ri(t  — u)dtdu 

+  Riit  —  dtdu  +  m(u)  Ri(t  —  u)  dtdu. 
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Performing  a  change  of  variable  by  substituting  y=zt-u,  the  following  is 
obtained; 

=  fo  -RiW  +  foM^)  Ri(y)  dydu 

+  /o'^Tn(u)  Ri^y)  dydu  +  J^m{u)  fj  e“^ Ri(y)  dydu 

=  fj  e"^'  Ri:,(t)  dt  +  fje~^'^m(u)  i?i(y)  dydu 

+  jJe-^'^m(u)fj_^e~^yRi(y)dydu  +  /■^e~^'^m(u)  fJ  Ri(y)  dydu 

=  fJ RUt)  dt  +  fje-^'^m(u)  fJ  R^  dydu 

+  /■^e"®“7n(u)  fJ  Riiv)  dydu 

^’(5)  =  fo  dt  +  ((/o°° e“*“m(u) du)  (/J i?i(y) dy))  .  (4.9) 

By  definition  the  following  holds: 

iT^*{s)  =  f^e~^^m(u)  du  . 

Define  the  partial  or  truncated  Laplace  transform  of  a  function,  f(t),  as 

ns,T)  =  ffe-“mdt.  (4.10) 

Thus,  the  following  holds; 

Rl{s,T)  =  So‘~“Mt)dt  . 
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Using  this  notation,  equation  (4.9)  becomes  the  following: 


A-(s)  =  Rl(s,  T)  +  (m'(s)  +  Rl(s.  T)) 

A-W  =  iJ2(s,T)(l  +  m-(s))  .  (4.11) 


Substituting  the  results  of  equation  (4.5)  into  (4.1 1)  the  following  is  obtained: 


Simplifying,  the  desired  result  is 

Equation  (4.12)  gives  the  Laplace  transform  of  the  availability  function  for  a 
component  age  replacement  preventive  maintenance  policy.  In  this  form,  it 
does  not  require  the  Laplace  transform  of  the  renewal  density.  It  depends  on 
the  partial  Laplace  transform  of  the  component's  survivor  function  and  the 
Laplace  transform  of  the  probability  density  function  for  the  time  between 
renewals  for  the  referenced  stochastic  process. 

It  is  assumed  that  the  distribution  function  for  the  component's  lifetime  is 
known  and  therefore  the  survivor  function,  Riit),  equivalent  to  (1-Fi(<))  is  also 
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known.  Also,  the  probability  density  functions  for  the  time  to  complete  a  repair 
due  to  component  failure  as  well  as  the  time  to  complete  a  preventive 
maintenance  action  are  known.  However,  the  probability  density  function  for 
the  time  between  renewals,  referred  to  as  h{t),  is  not  known.  A  renewal  theory 
approach  is  applied  to  develop  the  Laplace  transform  of  the  probability  density 
function  for  the  time  between  renewals. 


4.4  Laplace  Transform  for  the  Time  Between  Renewals 

The  notation  used  in  this  development  is  first  presented  in  Section  3.1 
and  for  convenience  is  reiterated  below.  The  stochastic  process  under 
consideration  is  depicted  by  Figure  3.1  in  Section  3.1. 

Notation: 

fi{t)  =  known  component  lifetime  density 
Flit)  =  known  component  lifetime  distribution,  jQfL(t)dt 
Riit)  =  known  component  survivor  function,  1  -  Fi(t) 
g^{t)  =  known  repair  time  density  after  component  failure 
5p(f)  =  known  preventive  maintenance  time  density 

T  =  age  replacement  preventive  maintenance  period,  i.e.  if  a 
component  survives  T  time  units  then  it  is  preventively 
maintained. 
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The  following  assumptions  are  made  concerning  this  stochastic  process: 

(1)  All  probability  densities  (lifetime,  failure  repair,  and  preventive 
maintenance  action)  are  known  or  may  be  estimated  from  available  data. 

(2)  A  single  component  is  modeled  and  is  considered  to  be 
continuously  demanded  when  available.  If  the  component  is  not  continuously 
demanded,  then  it  continues  to  age  according  to  the  lifetime  density  during 
periods  of  availability  whether  demanded  or  not. 

(3)  Repairs  as  a  result  of  component  failure  and  preventive 
maintenance  actions  are  perfect.  At  their  completion,  the  component  is  "as 
good  as  new"  and  thus  these  points  are  renewal  points  in  the  process.  At  the 
beginning  of  each  cycle,  the  component  exhibits  identical  probabilistic 
behavior  as  a  "brand  new"  component. 

(4)  Upon  completion  of  a  repair  due  to  component  failure  or  a 
preventive  maintenance  action,  the  component  is  immediately  placed  back  into 
operation. 

(5)  Upon  completion  of  a  repair  due  to  component  failure  or  a 
preventive  maintenance  action,  the  time  counter  for  age  replacement 
preventive  maintenance  policy  is  reset  to  zero. 

Two  types  of  events  or  cycles  occur  in  this  stochastic  process.  The  first 
is  termed  a  "failure  cycle".  In  this  event,  the  component  fails  before  time  T  and 
undergoes  a  stochastic  repair  period  subject  to  the  probability  density, 

The  second  event  is  termed  a  "preventive  maintenance  cycle".  In  this  event, 
the  component  survives  until  time  T  and  then  undergoes  a  stochastic 
preventive  maintenance  period  subject  to  the  probability  density,  To 
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develop  the  probability  density  function  for  the  time  between  renewals  in  this 
process,  an  analysis  of  each  cycle  is  first  performed.  The  Laplace  transform 
for  the  probability  density  function  of  the  time  between  renewals  within  each 
cycle  is  developed.  Then,  by  conditioning  on  the  time  of  component  failure,  an 
expression  for  the  Laplace  transform  for  the  time  between  renewals,  h(t),  is 
developed. 


4.4.1  Failure  Cycle 

In  this  event,  the  component  fails  before  time  T.  Thus,  if  t  is  the 
continuous  time  variable,  then  t  <J.  Recall  that  we  assume  that  the 
component  lifetime  density  function  is  known.  However,  in  this  cycle  the 
lifetime  density  function  is  truncated  at  time  T.  This  requires  the  component 
lifetime  density  function,  /^(f),  to  be  normalized  by  the  factor,  Fi(T).  In  fact, 

Fl(T)  is  the  probability  that  the  component  fails  before  time  T.  The  truncated 
lifetime  density  function  is  fort  <  T  and  Q otherwise.  Define  the 

probability  density  function  for  the  time  between  renewals  in  the  failure  cycle 
as  hf{t).  This  density  function,  hf{t),  is  the  convolution  of  the  truncated  lifetime 
density  function,  and  the  failure  repair  time  probability  density  function, 

g^(t).  This  is  expressed  mathematically  as 

A/W  =  (^)  ®  9.W  .  (4.13) 


58 


where  the  cx)nvoIution  operator  is  designated  by  the  symbol,  (g) .  Taking  the 
Laplace  transform  of  both  sides  of  equation  (4.13)  and  applying  the 
convolution  theorem  for  truncated  Laplace  transforms  yields  the  following: 


=  (/o’''""  ^ 

'•J  =  fdr)  (//' fdt)  dt)  . 

Note  that  the  Laplace  transform  of  the  truncated  lifetime  density  is  also 
truncated  at  T.  This  is  reflected  in  the  upper  limit  of  the  integral  since  the 
interval  for  t  is  [  0,T  ]  in  this  cycle.  The  component  fails  before  T  time  units 
have  elapsed. 


4.4.2  Preventive  Maintenance  Cycle 

In  this  event,  the  component  survives  until  time  T.  Thus,  if  t  is  the 
continuous  time  variable,  then  f  =  T.  In  this  cycle  the  operation  time  of  the 
component  is  a  constant  time  period  equal  to  T  time  units.  Probabilistically, 
this  constant  time  period  may  be  treated  as  a  Dirac  Delta  function  (or  Unit 
Impulse  function)  for  the  probability  density  function.  This  function  is 
designated,  6(t  -  a)  and  has  a  singularity  of  infinite  value  at  f  =  a  and  is  equal 
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to  zero  at  all  other  values  of  t.  The  Dirac  Delta  function  for  this  constant  time 
period,  6{t  -  T),  has  the  following  properties  ([102],  pp.  8-9): 

(1)  J^6it-T)dt=l 

(2)  J^6{t  -  T)  G(t)  dt  =  G(J)  for  any  continuous  function. 

Property  1  gives  the  Dirac  Delta  function  the  primary  characteristic  of  a 
probability  density  function.  Property  2  assists  in  taking  the  Laplace  transform 
of  the  Dirac  Delta  function.  The  Laplace  transform  for  the  Dirac  Delta  function 
for  this  constant  time  period  is 

£■>[«(( -T)l  =  /„”«(( -T)e-*‘d( 

Define  the  probability  density  function  for  the  time  between  renewals  in  the 
preventive  maintenance  cycle  as  hp{t).  This  density  function  is  the  convolution 
of  the  Dirac  Delta  function,  6{t  -  T),  and  the  preventive  maintenance  time 
probability  density  function,  gp{t).  This  is  expressed  mathematically  in  the 

following; 


hp{t)  =  (6(t  -  T))  (8>  Qpit), 


(4.14) 
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where  the  convolution  operator  is  designated  by  the  symbol.  <s> .  Taking  the 
Laplace  transform  of  both  sides  of  equation  (4.14)  and  applying  the 
convolution  theorem  for  Laplace  transforms  yields 


(4.15) 


4,4.3  Conditioning 

To  obtain  the  complete  expression  for  the  Laplace  transform  of  the 
probability  density  function  for  the  time  between  renewals,  h'is),  for  the 
stochastic  process  comprised  of  randomly  occurring  failure  cycles  and 
preventive  maintenance  cycles,  we  condition  on  the  length  of  time  the 
component  has  been  in  operation  at  the  time  of  component  failure  (refer  to 
Figure  3.1).  Specifically,  two  conditions  shall  be  used.  These  are: 

(1 )  Component  failure  before  the  component  has  been  in  operation  for 
T  time  units  since  the  last  service  action  (foper  <  T,  where  toper  is  elapsed  time 
of  component  operation  since  the  last  service  action).  A  service  action  is 
either  a  repair  due  to  component  failure  or  a  preventive  maintenance  action. 
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(2)  Component  survival  until  the  component  has  been  in  operation  T 
time  units  since  the  last  service  action  (ioper  >  T,  where  toper  is  elapsed  time  of 
component  operation  since  the  last  service  action). 

The  /i*(s)  may  be  found  by  using  the  following  conditional  probability 
statement; 


h  (s)  =  /i*(s)  Pro6{toper  <  T}  +  Pro6{toper  >  T] 

^ooer^'r  ^opar^'T 


Under  the  condition,  toper 
cycle  event.  Thus, /i*(s)| 


<  T,  the  stochastic  process  experiences  a  failure 
,  is  as  shown  in  equation  (4.13).  Under  the 

^oper'^T 


condition,  toper  >  T,  the  stochastic  process  experiences  a  preventive 
maintenance  cycle  event.  Thus,  h*{s)  ,  is  as  shown  in  equation  (4.14). 

^ODer^*r 


The  probability  that  the  component  fails  before  it  operates  T  time  units, 
Pro6{toper  <  T},  is  Pl(T).  The  probability  that  the  component  survives  at  least 
T  time  units,  Pro6{toper  >  T},  is  Rl{T)  =  1  -  Pi(T).  Substituting  these  into 
equation  (4.15),  the  following  expression  for  the  Laplace  transform  of  the  time 
between  renewals  for  the  referenced  stochastic  process  is 


h-(B)  =  !U^.  l)9'M+e-^RLmg;(s)  ,  (4.16) 


where  /2(s,  T)  is  the  partial  Laplace  transform  of  as  defined  by  equation 
(4.10). 
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4.5  Availability  Model  Conclusion 


The  availability  model  for  the  age  replacement  preventive  maintenance 
is  complete.  An  expression  is  given  for  the  Laplace  transform  of  the  availability 
function  in  terms  of  the  Laplace  transforms  of  the  component's  survivor 
function  and  the  probability  density  function  for  the  time  between  renewals. 
Restated,  ^’(s)  is 


A’(s)  = 


Rds.T) 

l-h\s)  ’ 


(4.17) 


where  h*{s)  is  shown  in  equation  (4.16). 


4.5.1  Model  Issues 


Three  main  issues  affect  the  use  of  this  model  to  obtain  an  expression 
for  A{t).  These  issues  include  finding  partial  (truncated)  and  normal  Laplace 
transforms,  finding  an  estimated  inverse  Laplace  transform  for  A*{s)  and  the 
numerical  integration  of  A{t)  to  find  the  value  of  T  that  maximizes  the  average 
availability  of  a  time  interval  [0,  t]. 
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4.5. 1.1  Partial  and  Normal  Laplace  Transforms 


The  first  issue  involves  finding  the  partial  Laplace  transforms  of  the 
lifetime  density,  /^(s,  T)  and  the  partial  Laplace  transform  of  the  component's 
survivor  function,  Rl{s,  T).  For  the  most  interesting  types  of  lifetime  densities 
(i.e.,  those  exhibiting  the  ability  to  model  increasing  hazard  or  "component 
aging"),  special  power  expansion  techniques  are  employed.  The  lifetime 
density  to  be  used  in  this  case  is  the  Weibull  density  function  defined  [95]  as 

flit)  -  ,  t  >  0, 

and  Flit)  =  1  -  6-“*', 

where  a  is  the  scale  parameter  and  p  is  the  shape  parameter.  For  values  of 
/3  >  1,  the  Weibull  lifetime  density  exhibits  strictly  increasing  hazard. 

Consider  the  Rayleigh  density  function  which  is  a  special  case  of  the 
Weibull  with  /5  =  2.  In  this  case,  the  hazard  function,  z{t),  is  defined  as 
follows; 

^  Rdt) 

Therefore,  for  the  Rayleigh  lifetime  density  function,  the  hazard  function  is 
strictly  increasing  in  a  linear  fashion.  Practically,  this  family  of  distributions 
models  aging  components  but  does  not  model  a  "wear  out"  period  very  well 
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since  the  hazard  function  increases  at  same  rate  over  all  time.  Wear  out 
periods  may  be  modeled  using  Weibull  lifetime  densities  with  values  of  /3  >  2. 
In  these  cases,  the  force  of  component  mortality  increases  exponentially  with 
time.  Analysis  here  proceeds  with  the  Rayleigh  distribution  and  the  application 
of  power  series  expansions  to  find  fl{s,  T)  and  RKs,  T).  These  expansions 
are  truncated  to  provide  estimates  for  the  preceding  terms. 

In  the  event  that  the  Weibull  distribution  is  used  to  model  the  repair  time 
and  preventive  maintenance  time  distributions,  a  similar  approach  is  followed. 
The  difference  in  this  case  is  that  these  service  times  are  not  truncated.  Thus 
we  must  find  an  estimate  for  the  full  Laplace  transform.  Power  series 
expansions  may  be  applied  to  find  g^{t)  and  gp{t).  If  successful,  these 

expansion  are  truncated  to  provide  estimates  for  the  preceding  terms. 


4.5. 1.2  Laplace  Transform  Inversion  Estimation 

The  complexity  of  the  expression  obtained  for  A* (s)  will  increase  as 
more  realistic  distributions  are  used  for  the  component's  lifetime  density  as 
well  as  the  repair  and  preventive  maintenance  action  times.  Due  to  this 
complexity,  it  is  unlikely  that  an  exact  inversion  of  A*(s)  may  be  accomplished. 
Therefore  methods  to  estimate  the  time  domain  function,  A{t),  are 
investigated.  One  family  of  methods  provides  an  estimate  of  the  value  of  the 
time  domain  function,  A{t),  at  specific  points  in  time  from  the  A*(s).  One 
approach  would  be  to  iteratively  run  a  numerical  inversion  technique  over 
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many  values  of  time  and  then  fit  a  function  to  the  values  to  obtain  an  estimate 
of  A{t).  Unfortunately,  this  methodology  would  become  very  time  consuming 
and  cumbersome.  However,  Cheng  et  al  [66]  have  implemented  several  of  the 
well  known  Laplace  transform  numerical  inversion  techniques  in  Mathematica. 
Their  focus  is  on  the  numerical  solution  of  linear  partial  differential  equations 
relating  to  groundwater  and  reservoir  engineering  applications.  Mathematica 
is  used  to  compile  numerical  results  for  functions  with  known  Laplace  inverses 
for  comparison  purposes.  Due  to  the  analytical  capabilities  of  Mathematica, 
this  analysis  uses  these  inversion  estimation  techniques  in  a  somewhat 
different  manner.  In  the  analysis  presented  here,  the  time  parameter,  t,  is  left 
as  a  parameter  in  the  inversion  techniques.  The  result  from  Mathematica  is  an 
estimate  of  the  time  domain  function,  A{t),  which  is  an  extension  from  a  single 
estimated  value  of  the  time  domain  function  at  a  specific  time,  t. 


4.5.1 .3  Numerical  Integration  Search  to  Find  Values  of  T  to 
Maximize  Average  Availability 

After  obtaining  an  estimate  of  A{t),  a  numerical  integration  search  is 
attempted  for  a  specific  value  of  t,  to  find  the  value  of  T  that  maximizes 
average  availability  over  the  time  interval  [0,  t].  This  numerical  search  is 
performed  in  Mathematica.  The  results  are  compared  to  the  values  of  T  for  the 
infinite  time  horizon  case. 
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Additionally,  methods  for  obtaining  a  direct  functional  form  of  the 
average  availability  with  respect  to  time  are  investigated. 
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CHAPTER  5 


THE  EXPONENTIAL  CASE 


5,1  Exponential  Case  Description 


In  this  case,  the  component's  lifetime  as  well  as  the  time  required  to 
perform  repair  upon  component  failure  and  the  time  required  to  perform 
preventive  maintenance  on  the  component  are  all  considered  to  have 
exponential  probability  density  functions  (pdf).  These  exponential  density 
functions  are  assumed  to  be  unique  and  independent  of  each  other.  These 


functions  are  defined  as 

II 

> 

I 

> 

(the  component  lifetime  pdf). 

= 

(the  repair  time  pdf),  and 

9pii)  =  y-p  e-V 

(the  preventive  maintenance  time  pdf) 

This  particular  case  is  selected  to  provide  validation  examples  for  the 
proposed  model.  The  hazard  function,  z{t),  is  defined  as  follows  [7]: 


For  the  exponential  lifetime  density,  z(t)  is  equal  to  A,  a  constant.  Thus,  as  the 
component  ages,  the  hazard  remains  constant.  Therefore  the  "memoryless" 
property  of  the  exponential  pdf  applies.  Given  that  the  component  survives  t 
time  units,  the  probability  that  the  component  fails  in  the  next  instant  of  time 
remains  constant  for  all  t.  In  essence,  the  component  does  not  exhibit  aging 
since  whether  it  has  survived  one  time  unit  or  one-hundred  time  units,  the 
probability  of  it  failing  in  the  next  instant  of  time  is  the  same.  Thus,  given  the 
component  has  survived  t  time  units,  it  is  always  as  good  as  new.  Barlow  and 
Proschan  [7]  term  this  constant  hazard. 

Logically,  it  is  easy  to  argue  that  if  the  component's  lifetime  is  modeled 
in  this  manner  then  preventive  maintenance  should  never  be  performed. 

There  would  not  exist  a  component  age  at  which  one  would  desire  to  purposely 
remove  the  component  from  operation  to  perform  preventive  maintenance. 
Recall  that  given  the  component  has  survived  t  time  units,  it  is  as  good  as  new. 
Thus  for  the  exponential  case,  the  replacement  period  maximizing  average 
availability  should  be  infinity  or  in  other  words,  no  preventive  maintenance 
should  be  performed.  In  this  case,  the  process  simplifies  into  an  alternating 
renewal  process  where  the  component  is  repaired  upon  failure. 

The  model  developed  in  Chapter  4  is  used  to  analyze  the  exponential 
case  described  above.  This  case  offers  a  method  to  assist  in  validating  the 
model  since  the  correct  outcome  is  known. 
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5.2  Development  of  the  Availability  Function  for  the  Exponential 
Case 


Recalling  equations  (4.16)  and  (4.17)  from  Sections  4.4.3  and  4.5 
respectively,  the  Laplace  transform  of  the  availability  function.  >l*(s),  is 


where  h-(s)  =  fl(s,  T)s;(s)  +  e->'JZi(T)s;(s) 

and  /■(s,T)  =  /„V'V(f)<i(  . 


Applying  these  equations  to  the  referenced  exponential  case,  the  following  is 
obtained: 


dt 


A*(e) 


±  (.-TWA))) 


A*(e)  = 


_ if^p+s)  {fJr+s)  -  (/Jp+s)  (lJ.r+s)  _ 

(A+s)  i/ir+s)  iHp+s)-Xn^{Hp+s)+(^Xfji^ifij,+s)-HpiX+s)  (/ir+s)) 


(5.1) 
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5.3  Laplace  Transform  Inversion  to  Obtain  Ait^ 


Inversion  of  this  function  shown  in  equation  (5.1 )  to  the  time  domain  is 
not  simple.  One  simplifying  condition  is  to  note  that  the  numerator  has  two 
terms  that  only  differ  by  a  factor  of  This  allows  the  inversion  to  be 

broken  down  into  the  inversion  of  two  terms.  This  factor  e“^^®+''^may  be  stated 
as  in  which,  e'^^,  is  a  constant.  A  useful  property  of  inverting  Laplace 

transforms  (see  [98])  is 

£-i[g-T5g-TA^*(^)]  ^  g-TA  (unitstep[f-T])  F(f-T)  ,  (5.2) 

where  F{u)  =  £'^[/(s)] 

and  Unitstep[t-Tl  =  {?;^;:J^„»  . 

This  allows  the  inversion  to  proceed  with  only  the  first  term  and  then  use  the 
rules  shown  above  to  obtain  the  inverse  of  the  second  term.  We  then  sum  the 
results  to  obtain  the  A{t). 

The  work  now  lies  with  inverting  the  function  shown  below. 


g-T(5+A) 


1  (6.3) 
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Using  Mathematica  (the  Apart  function)  this  term  may  be  broken  down  into  the 
sum  of  the  following  fractions: 


P-1  (  (^r+g)  , 

\  s{fi^+\+s)  j 

£-1/  _ (Mr+s)(^^/Jp-/Z^A+/ipA+/ip5) _ 

\  s(A‘r+A+s)(/irA-/irMp-A'pA-MpS+(e(^+'')^(;ir/ip+/ipA+(/i^+/ip+A)s+s2) 


(5.4) 


The  first  term  in  this  sum  may  be  inverted  exactly  and  is 

£-1  (  ..itr+s)  \  p-t(/X,+A) 

\  s(Mrr+A+s)  )  Mr+A  /^r+A 

The  time  domain  function  for  the  first  term  in  equation  (5.4)  is  identical  to  the 
time  domain  function  for  the  alternating  renewal  process  where  repair  upon 
failure  is  executed  (i.e.  preventive  maintenance  is  not  performed).  The 
problem  now  is  to  invert  the  second  term  of  the  summation  shown  in  equation 
(5.4).  Note  that  this  term  in  equation  (5.4),  which  is  designated  as  g^is),  may 
be  represented  as 

=  S)  ■ 

Since  q{s)  is  a  second  order  polynomial  with  respect  to  s.  and  of  the  form 
Cl  +  C2S  +  C3s^,  we  may  apply  standard  rules  for  inverting  the  Laplace 
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transform.  If  G2{t)  represents  the  £‘*{^2(3)}  and  P{t)  represents  the  E'HirT} 

p(5j  ^ 

[97],  we  obtain 


Gj(()  =  CiP(t)  +  c,  |[P(t)]  +  C3  ^[P(<)] .  (5.6) 


This  decomposition  process  leaves  the  task  of  inverting  However,  this 

term  is  successfully  inverted  exactly.  An  estimate  of  the  inverted  function  was 
obtained  using  a  numerical  inversion  package  obtained  from  Wolfram 
Research,  Inc.,  distributors  of  Mathematica. 


5.3.1  Laplace  Transform  Inversion  Estimation 

A  numerical  inversion  package  for  Mathematica  written  by  Cheng,  et  al, 

[66]  is  implemented.  The  Stehfest  method  [67  and  98]  is  used  to  invert  the 
( ^ )  term  referenced  in  the  previous  section.  The  Stehfest  method  for 

approximating  inverse  Laplace  transforms  is  designed  to  provide  an  estimate 
of  the  time  domain  function  at  a  specific  value  of  time  from  the  Laplace 
transform  function.  For  this  reason,  the  technique  is  termed  a  "numerical" 
inversion  package.  The  power  of  Mathematica  allows  the  time  specification  to 
be  left  as  a  parameter.  This  allows  Mathematica  to  analytically  solve  for  the 
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time  domain  function  and  hence  give  us  an  estimate  of  the  time  domain 
function,  A{t). 

The  Stehfest  method  [67  and  98]  obtains  an  estimate  of  the  time  domain 
function  by  defining  the  time  domain  function,  F{t),  as 

n=l 


where 


£[F(f)]  =  f*{u)  and 


Cn  — 


Tnin{n,  y) 


kr(2k)\ 


( f  -A:)!  k\  (A:-l)!  (n-fc)!  (2fc-n)! 


Tests  were  performed  within  Mathematica  to  verify  the  validity  of  the 
Stehfest  method.  The  tests  consisted  of  performing  the  numerical  inversion  for 
a  specific  value  of  time  for  a  function,  f*{s).  The  result  of  this  test  provided 
time  domain  values  for  F{t)  (note  that  £'^[/*(s)]  =  F(t)) .  The  next  test 
consisted  of  performing  the  numerical  inversion  with  the  time  specification  left 
as  the  parameter,  t.  The  result  was  an  estimate  of  the  function,  F(t)  (as  a 
function  of  time).  The  same  values  used  in  the  previous  test  were  evaluated 
by  the  estimated  time  domain  function  and  were  found  to  be  the  same  as  those 
obtained  by  direct  numerical  inversion.  The  test  was  successful  and  the 
"numerical"  inversion  package  was  implemented  to  provide  an  estimate  of  the 
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time  domain  function.  The  results  of  these  tests  appear  in  Appendix  1 .  This 
appendix  provides  the  annotated  Mathematica  input  statements  and  output 
used  in  obtaining  these  results. 

Using  the  relationship  shown  in  equation  (5.2),  an  estimate  of  the 
inversion  of  the  first  term  in  equation  (5.1 )  shown  in  equation  (5.2)  is  obtained. 
Finally  an  estimate  of  the  availability  function,  A(t),  is  obtained  by  applying 
equation  (5.2).  Appendix  2  shows  the  Mathematica  statements  used  to  find 
this  estimate  of  A{t)  for  the  exponential  case.  This  estimated  function  is  long 
and  complicated  and  is  left  as  a  function  in  Mathematica.  If  desired,  the 
estimated  function  obtained  from  the  Stehfest  method  may  be  observed  within 
Mathematica.  Note  that  in  these  results,  T,  the  age  replacement  period  is 
included  as  a  parameter.  Thus  the  estimated  availability  function  obtained, 
A(t,  T),  is  actually  a  function  of  both  continuous  time,  t  and  the  age 
replacement  period,  T. 


5.4  Model  Validation 


Two  methods  are  used  to  validate  the  estimated  availability  function, 
A(t),  obtained  as  detailed  in  Section  5.3.1.  The  first  method  involves 

evaluating  the  value  of  the  availability  function  at  very  small  values  of  time. 
Since  the  estimated  function  contains  factors  of  (|),  evaluating  f  =  0  is  not 
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possible.  For  the  case  where  A=1.,  p  =  1000.  and  5  =  5.,  the  estimated 
availability  function  is  evaluated  for  t  =  0.001  for  various  values  of  T,  the  age 
replacement  policy.  In  all  cases  the  results  are  very  close  to  unity.  For 
example,  with  T  =  1.,  the  function  evaluated  at  f  =  0.001  was  0.999363.  Also 
for  T=10.,  the  function  evaluated  aXt  =  0.00000001  was  0.99999999.  These 
values  are  sufficiently  close  to  unity  and  provide  availability  model  validation 
for  small  values  of  t. 

The  second  technique  used  for  validating  the  model  is  evaluation  of  the 
availability  for  large  values  of  T  at  large  values  of  continuous  time.  t.  In  this 
case,  as  T-kx),  the  stochastic  process  becomes  a  repair  upon  component 
failure,  i.e.  no  preventive  maintenance  is  accomplished.  For  evaluating  this 
case,  we  set  A=1 .,  p  =  1000.,  6  =  b.  and  T  =  1000.  The  estimated  availability 
function  is  evaluated  at  t  =  500.  The  value  obtained  \sA(t  =  500.)  =  0.999001. 
The  limiting  availability  for  this  case  is  given  by 


A  =  Hm  A(t)  , 

t—*oo 

and  for  the  repair  upon  failure  case  (alternating  renewal  process)  is  equivalent 
to 

_  E[lifetime] 

E[lifetime]-{'E[repairtime] 

=  0-999001  ■ 
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Thus  for  the  values  shown  above  the  model  estimates  the  limiting  availability 
for  the  repair  upon  failure  only  process  exactly. 

These  validation  methods  basically  show  the  validity  of  the  model  at  the 
end  points.  The  first  technique  evaluated  the  model  as  while  the  second 
technique  considered  the  case  where  both  t  and  T  grew  large.  A  major 
implication  of  this  result  is  the  development  of  an  expression  for  the  availability 
function  for  the  component  in  this  process.  This  availability  function,  A(t,  T),  is 
a  function  of  both  time,  t  and  the  replacement  period,  T.  The  literature  search 
did  not  show  this  result  in  past  research.  Thus  a  method  to  evaluate  the 
model's  validity  between  the  extremes  discussed  above  has  not  been  devised. 
However,  the  validation  results  presented  are  encouraging  and  analysis  of  the 
availability  model  will  proceed. 


5.5  Maximizing  Average  Availability  Over  Finite  Time  Intervals 


The  average  availability,  Aavgir),  is  defined  as  [7] 


Aavgir)  =  A{t)dt  . 
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The  desired  result  is  to  find  a  value  for  T,  the  age  replacement  policy,  that 
maximizes  the  average  availability  over  a  finite  time  horizon,  t.  Two  methods 
to  find  the  value  of  T  are  detailed  in  Chapter  3  .  These  methods  were  an 
integral  equation  search  method  and  using  Laplace  transform  inversion  rules 
to  obtain  an  analytical  estimate  of  Aavg{r)  directly.  The  latter  method  is 
implemented  and  compared  to  cases  in  which  the  first  method  is  used. 

Function  values  correspond  well  and  therefore  the  method  of  using  Laplace 
transform  rules  to  obtain  an  analytical  estimate  of  Aavgir)  is  not  used. 

Mathematica  is  used  to  perform  searches  for  T  that  maximizes  the  value 
of  Aavgir)  for  specific,  finite  values  of  r.  An  iterative,  "double  Do  loop"  is 
implemented  through  the  Table  function  within  Mathematica.  The  algorithm 
iterates  through  values  of  finite  time  horizon,  r.  At  each  value  of  r,  Aavgir) 
function  values  are  calculated  for  several  values  of  T.  Meaningful  Aavgir) 
function  values  only  exist  for  the  case  where  T  <  r.  The  finite  time  horizon,  r, 
is  the  operational  period  of  the  component's  lifecycle.  Thus  if  T  >  r,  then  no 
preventive  maintenance  is  performed.  This  analysis  is  constrained  to 
evaluating  values  of  Aavgir)  for  which  T  <  r.  If  the  maximum  value  of  Aavgir) 
occurs  at  T  =  T,  then  the  age  replacement  period  is  equivalent  to  the  finite  time 
horizon.  Therefore,  no  preventive  maintenance  would  be  performed. 

For  the  exponential  case,  as  explained  earlier,  no  preventive 
maintenance  should  be  performed.  This  is  true  regardless  of  the  time  horizon. 
An  additional  concern  is  the  relationship  between  the  mean  repair  time  due  to 
failure  and  the  mean  preventive  maintenance  time.  If  the  mean  preventive 
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maintenance  time  exceeds  the  mean  repair  time,  Barlow  and  Hunter  [12],  have 
shown  that  no  preventive  maintenance  should  be  performed.  The  objective  in 
obtaining  this  result  is  to  maximize  the  limiting  availability.  The  ratio,  6,  defined 
to  be  the  mean  repair  time  divided  by  the  mean  preventive  maintenance  time  is 
always  be  greater  than  one  {6  >  1).  Furthermore,  the  ratio,  p,  defined  to  be  the 
mean  lifetime  divided  by  the  mean  repair  time  is  characteristically  large. 

Values  of  10  to  10^  may  be  typical  depending  on  the  class  of  the  component 
(i.e.  mechanical,  electronic).  For  demonstration  purposes,  smaller  values  of  p 
are  used  to  exaggerate  values  of  availability  for  ease  of  analysis. 


5.5.1  Exponential  Model  Results 


The  purpose  of  these  results  is  to  provide  additional  model  verification. 
The  correct  result  for  all  values  of  r,  is  to  obtain  a  model  result  that  shows  the 
maximum  Aavg{r)  to  be  attained  when  preventive  maintenance  is  not 
performed.  It  is  expected  that  this  result  should  be  evident  by  Aavg{r)  being 
maximized  when  T  =  r. 

The  following  exponential  model  is  evaluated.  The  mean  lifetime, 
i  =  1.;  p  =  2.,  thus  the  mean  repair  time,  —  =  0.5  .  The  ratio  of  mean  repair 

A  /if 

time  to  mean  preventive  maintenance  time  is  5  =  5.  Thus,  the  mean 
preventive  maintenance  time,  —  =  0.1  .  The  results  are  shown  below  for 
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values  of  r  between  0.5  and  5.  Values  of  T  include  the  upper  end  point.  The 
results  are  shown  in  Table  5.1 .  Observe  that  the  maximum  value  of  Aavgir) 
occurs  at  the  point  where  T  =  r.  In  some  instances,  the  number  of  digits 
shown  does  not  reflect  this  maximization.  However,  using  sixteen  digit 
accuracy  within  Mathematica  reflected  Aavgij)  being  maximized  at  T  =  r.  This 
result  confirms  that  for  the  exponential  case,  preventive  maintenance  should 
not  be  performed  to  maximize  Aavgir)  for  finite  time  horizons.  This  provides 
additional  evidence  supporting  the  validity  of  the  model. 
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CHAPTER  6 


THE  WEIBULL  -  EXPONENTIAL  CASE 


6.1  Introduction 


In  this  case,  I  examine  the  age  replacement  preventive  maintenance 
stochastic  process  with  a  Weibull  component  lifetime  probability  density 
function  and  exponential  repair  and  preventive  maintenance  (p.  m.)  time 
probability  density  functions.  These  density  functions  are  defined  as 

f^(t)  =a(3 (component  lifetime  pdf), 

g^{t)  =  (repair  time  pdf), 

g^(t)  =  Up  (preventive  maintenance  time  pdf),  and 

T  (age  at  which  component  receives 

preventive  maintenance). 

The  stochastic  process  consists  of  randomly  occurring  failure  cycles  and 
preventive  maintenance  cycles.  A  failure  cycle  occurs  when  the  component 
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fails  before  it  ages  T  time  units.  The  Weibull  distribution  governs  the  random 
failure  of  the  component.  Given  that  the  component  fails  before  time  T  (before 
T  time  units  have  elapsed  since  the  last  service  action),  the  component 
immediately  begins  failure  repair.  This  random  repair  period  is  exponentially 
distributed.  If  the  component  survives  T  time  units  (since  the  last  repair  or 
preventive  maintenance  action)  then  the  component  is  immediately  removed 
from  service  and  begins  preventive  maintenance  service.  The  time  required 
for  preventive  maintenance  is  random  and  is  exponentially  distributed.  The 
exponential  distributions  modeling  the  preventive  maintenance  action  and 
repair  times  are  unique  and  independent  of  each  other. 


6.2  Application  of  Model 


Recalling  equations  (4.16)  and  (4.17)  from  Sections  4.4.3  and  4.5,  the 
Laplace  transform  of  the  availability  function,  A*(s),  is 


A*(5) 


i?l(g.T) 

l-h\s)  ’ 


where  h’is)  =  fl(s,  T)  g-{s)  +  g;(s) 


and 


/I(5.  T)  = /oV*ViW  . 
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Applying  these  equations  to  the  referenced  Weibull-exponential  case  requires 
the  development  of  expressions  for  the  partial  Laplace  transform  for  the 
Weibull  probability  density  function,  fl(s,  T),  as  well  as  for  the  Weibull  survivor 
function,  R^is,  T).  The  Laplace  transforms  for  the  exponentially  distributed 
repair  and  preventive  maintenance  times  are  not  truncated  since  the  range  on 
time  variable  in  these  cases  is  0  <  f  <  oo.  These  transforms  are  easily 
evaluated. 


6.2.1  Laplace  Transform  of  the  Weibull  Probability  Density  Function 

In  the  form  given  in  Section  6.1 ,  the  Laplace  transform  of  the  Weibull 
probability  density  function  may  not  be  found.  However,  the  Weibull  pdf  may 
be  expanded  into  a  power  series.  The  Laplace  transform  of  the  power  series 
representation  of  the  Weibull  pdf  may  be  found.  Thus,  an  exact  Laplace 
transform  of  the  Weibull  pdf  is  generated  in  the  form  of  an  infinite  series. 

The  methodology  used  to  obtain  the  Laplace  transform  for  the  Weibull 
pdf  is  also  used  to  find  the  partial  Laplace  transforms  of  the  Weibull  pdf  and 
survivor  function. 
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6.2.1. 1  Power  Series  Expansion  of  the  Weibull  Probability  Density 
Function 


Recall  that  the  Weibull  probability  density  function  is  defined  as  [95] 

=  (6.1) 

Using  a  similar  methodology  as  that  shown  by  Lomnicki  [99],  this  function  may 
be  expanded  into  a  power  series.  First,  let  x  =  t^  and  expand  the  exponential 
term  of  the  function  as 


00 

n=0 

Substituting,  f^for  x,  equation  (6.2)  is  obtained. 

e-“‘'’  =  i:(-ir(^)  .  (6.2) 

n=0 


Substituting  this  MacLaurin  series  expansion  into  equation  (6.1),  the  power 
series  expansion  for  the  Weibull  probability  density  function  is 


OO 

SM  =  E(-ir‘ 

n=1 


(n-1)! 


(6.3) 
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6.2.1 .2  Laplace  Transform  of  the  Weibull  Power  Series  Expansion 


By  definition,  the  Laplace  transform  of  a  function,  £[/(f)]  is 


mt)]  = 


Applying  this  to  equation  (6.3),  the  following  is  derived: 


fli^)  =  /o°° 


m^)  =  ( E(-ir*  )  dt 

'  n=1  ' 


(6.4) 


Mathematica  is  used  to  solve  the  integral  in  equation  (6.4).  The  resulting 
infinite  series  representation  of  the  Weibull  density  function  Laplace  transform 
is 


n=l 


(6.5) 


where  r[z]  is  the  Gamma  function  (Abramowitz  &  Stegun  [100],  p.  255  and 
Mathematica  Reference  manual,  [96],  p.  365)  and  is  defined  as 

r[z]  =  .  (6.6) 
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The  result  shown  in  equation  (6.5)  is  also  useful  for  analyzing  the  age 
replacement  preventive  maintenance  stochastic  process  with  Weibull  repair 
and  preventive  maintenance  time  distributions.  Note  that  these  time 
distributions  are  not  truncated  and  may  range,  0  <  f  <  oo.  However,  the 
immediate  requirements  for  this  analysis  are  the  partial  Laplace  transforms  of 
the  Weibull  lifetime  probability  density  function  and  the  Weibull  survivor 
function.  Sections  6.2.2  and  6.2.3  develop  these  partial  transforms. 


6.2.2  Partial  Laplace  Transform  of  the  Weibull  Lifetime  Density 
Function 


In  this  case,  we  are  concerned  with  the  time  interval  for  component 
failure  within  the  failure  cycle.  This  time  interval  is  random  and  is  distributed 
according  to  a  Weibull  probability  density  function.  This  known  probability 
function,  as  shown  in  equation  (6.1),  has  range  0<t  <  oo.  However,  this  time 
period  is  always  less  than  T  time  units  according  to  our  age-replacement 
preventive  maintenance  policy.  This  fact  is  taken  into  consideration  when  the 
expression  for  .4*(s)  is  developed  for  the  general  case  in  Chapter  4.  This  time 
truncation  causes  the  Laplace  transform  to  be  truncated  and  is  referred  to  as  a 
partial  Laplace  transform.  Recall  the  partial  Laplace  transform  we  are 
concerned  with  in  general  is 


rds,T,=j;e-‘mdt. 


(6.7) 
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Substituting  the  power  series  representation  of  the  Weibull  probability  density 
function  shown  in  equation  (6.3)  into  equation  (6.7)  gives 


—  /  00 

^  n=1 


(n-1)!  J 


(6.8) 


Mathematica  is  used  to  solve  the  integral  in  equation  (6.8).  The  resulting 
infinite  series  representation  of  the  Weibull  density  function  partial  Laplace 
transform  is 


T\  y[nl3,sT] 

J L\  ’  /  v  /  (n— 


n=l 


(6.9) 


where  7[a,  2]  is  the  alternative  incomplete  gamma  function  (Abramowitz  and 
Stegun  [100],  p.  260  and  by  Wolfram  [96],  p.  365)  and  is  defined  as 


7(0,2]  =  Jo  e  ^ dt  (6.10) 

Equation  (6.10)  gives  an  exact  (  as  n-»oo)  series  representation  of  the  partial 
Laplace  transform  of  the  Weibull  lifetime  probability  density  function.  Practical 
use  of  this  infinite  series  requires  truncating  the  infinite  series  to  obtain  a  finite 
representation  of  the  Weibull  pdf  partial  Laplace  transform.  This  development 
provides  the  basis  for  obtaining  the  partial  Laplace  transform  for  the  Weibull 
survivor  function,  Rl{s,  T).  This  development  is  presented  in  Section  6.2.3. 
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6.2.3  Partial  Laplace  Transform  of  the  Weibull  Survivor  Function 


The  survivor  function  for  a  component's  lifetime,  Rl{u),  is  defined  as 
the  probability  that  the  component  survives  at  least  a  period  of  time  u.  This  is 
stated  using  probability  notation,  where  u  represents  the  length  of  the 
component's  life  and  Y  is  some  specified  period  of  time,  as 

=  PTob{t  >  r}  =  1  -  F{Y) . 

In  this  case,  the  Weibull  survivor  function  is  defined  by  [7]: 


Rdt)  =  e-“ 


This  survivor  function  may  be  expanded  into  an  equivalent  power  series  by  the 
same  method  employed  in  Section  6.2.1. 1  .  The  power  series  representation 
of  the  Weibull  survivor  function  is 


00  „ 

«!(<)= (6.11) 

n=0 

The  partial  Laplace  transform  of  the  survivor  function,  Ri{t),  is 

Rl(s,T)  =  fJe-‘RUt)dt.  (6.12) 
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Substituting  the  power  series  representation  of  the  Weibull  survivor  function 
shown  in  equation  (6.1 1)  into  equation  (6.12)  gives 

RU^,  T)  -  (  E(-ir  ^  ■  (6.13) 

^  n=0  ■  ^ 


Mathematica  is  used  to  solve  the  integral  in  (6.13).  The  resulting  infinite  series 
representation  of  the  Weibull  survivor  function  partial  Laplace  transform  is 


RIU.  T)  =  1 


OO 

E(-i) 


7[(l+n/3),sT] 

(„)!s(n/3+l) 


(6.14) 


Equation  (6.14)  gives  an  exact  (  as  n-»oo)  series  representation  of  the  partial 
Laplace  transform  of  the  Weibull  lifetime  survivor  function.  Practical  use  of  this 
infinite  series  requires  truncating  the  infinite  series  to  obtain  a  finite 
representation  of  the  Weibull  survivor  function  partial  Laplace  transform. 


6.3  Numerical  Results  for  the  Weibull  -  Exponential  Case 

The  availability  model  is  constructed  in  Mathematica.  The  infinite  series 
for  the  partial  Laplace  transforms  of  the  Weibull  lifetime  density,  /2(s,  T),  and 
the  Weibull  survivor  function  are  truncated  at  32  terms  each.  Then  the 
availability  expression  developed  in  Chapter  4  is  implemented  in  Mathematica 
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to  obtain  an  estimate  of  the  availability  function  Laplace  transform,  -4*(s).  This 
transform  is  inverted  in  Mathematica  using  the  numerical  inversion  package 
developed  by  Cheng,  et  al  [66],  The  Stehfest  method  was  chosen  and 
implemented  forn  =  6  and  as  described  earlier  is  used  to  find  an  analytical 
estimate  of  A{t,  T). 

A  specific  numerical  case  is  used  in  this  analysis.  The  lifetime,  repair 
and  preventive  maintenance  densities'  parameters  are  shown  below.  The 
parameters  for  the  repair  and  preventive  maintenance  densities  reflect  values 
of  p  =  2.  and  6  =  5. 

Weibull  lifetime  density;  a=  /3  =  2,  meanlife  =  12.5, 

characteristic  life  =  14.14 

Exponential  repair  time  density:  =  0.16,  mean  =  6.25 

Exponential  preventive  maintenance  density:  pp  =  0.8,  mean  =  1.25 


6.3.1  Model  Validation  for  the  Weibull  -  Exponential  Case 

Three  methods  are  employed  to  validate  the  model  for  A{t).  First,  the 
estimated  value  of  A{t)  near  f  =  0  is  evaluated  for  several  different  values  of  T, 
the  age-replacement  time  interval.  Since  the  model  assumes  the  component  is 
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working  at  f  =  0,  the  value  of  A{t)  for  small  values  of  t  should  be  approximately 
equal  to  one  regardless  of  the  value  of  the  age  replacement  policy.  T  (T>0). 

The  second  method  is  the  computation  of  an  estimate  for  the  limiting 
availability  value  for  the  case  in  which  T  is  set  to  the  optimal  age  replacement 
period  for  an  infinite  time  horizon,  T^.  Barlow  and  Hunter  [12]  develop  the 
analytical  equations  to  solve  for  and  the  limiting.availability.  A.  This 
analysis  is  carried  out  by  first  setting  T=T^  and  then  evaluating  the  average 
availability,  Aavg{t),  for  large  values  of  time.  These  values  should  converge  to 
the  value  of  A  obtained  analytically. 

In  the  third  method,  the  estimated  availability  function  is  used  to 
compute  values  of  average  availability  for  several  values  of  T,  near  and 
including  T^  for  large  values  of  continuous  time,  t  For  large  values  of  t,  we 
expect  that  the  estimated  availability  function  should  produce  an  average 
availability  local  maxima  for  T=  T^. 


6.3. 1.1  Model  Validation  Results  for  f  0 


These  values  are  evaluated  in  the  estimated  availability  function,  A(t, 

T),  using  Mathematica,  for  t  =  0.000001.  This  evaluation  is  performed  for 
values  of  T  between  0.1  and  750.  The  expected  results  are  values  of 
A(t  =  0.000001,  T)  close  to  unity  for  all  values  of  T.  Table  6.1  shows  the  results 
of  this  analysis.  In  each  case,  Mathematica  returned  a  value  of  1 .  According 
to  the  Mathematica  documentation,  the  default  accuracy  is  at  least  sixteen 
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Table  6.1  A(t=0.0G0001)  for  0.1  <T  <750 


T  (age  replacement  period) 

A(t=0. 000001) 

0.10 

1. 

0.25 

1, 

0.50 

1. 

0.75 

1. 

1.00 

1. 

2.00 

1. 

5.00 

1. 

10.00 

1. 

20.00 

1. 

50.00 

1. 

100.00 

1. 

200.00 

1. 

500.00 

1. 

750.00 

1. 
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decimal  places.  These  results  support  the  validity  of  the  model  for  small 
values  of  time. 


6.3. 1.2  Validation  Results  for  the  Infinite  Time  Horizon 


The  results  of  this  validation  test  are  compared  to  the  limiting 
availability,  A,  obtained  from  the  infinite  time  horizon  case.  This  value  of  A 
depends  on  the  optimal  age  replacement  period.  T^.  Barlow  and  Hunter  [12] 
develop  these  analytical  equations.  T^  is  found  by  solving  the  following 
integral  equation  for  T^  (see  [12]): 

(^(^oo)  lo’^  -  Fl{T^)  -  E[Repairtime]-E[PM.  time]  ’ 

where  z{t)  is  the  hazard  function  and  £;[  •  ]  is  the  expected  value  operator. 

The  limiting  availability,  A,  is  found  by  solving  the  following  equation  (see  [12]): 


^  =  7 - r - 7—  ,  VO.  lOJ 

1  +  lE[Repair  time] -E[P.M.  time])  z(T^) 

These  equations  apply  when  the  lifetime  density  function,  has  a  strictly 
increasing  failure  rate  >  0)  and  the  mean  time  required  to  repair  the 

component  is  greater  than  the  mean  time  required  to  perform  preventive 
maintenance  {E{RepairtiTne\  >  E[P .M .time\)[\2\. 
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The  values  defined  above  for  the  component  lifetime  density  as  well  as 
the  repair  time  and  preventive  maintenance  time  densities  are  evaluated  in 
equation  (6.15)  to  find  T^.  Mathematica  solves  this  equation  quickly  and 
produces  =  7.22176.  This  value  of  T^  is  used  in  equation  (6.16)  to  obtain 
A  =  0.734706. 

The  value  of  T^,  =  7.22176  is  used  in  the  estimated  availability  function 
to  find  values  of  the  average  availability  ,  Aavg(t),  for  large  values  of  t.  Barlow 
and  Proschan  [7]  show  that  the  limiting  average  availability  is  equivalent  to  the 
limiting  availability.  Mathematically,  this  is 

lim  Aavg{t)  =  lim  -fjA(t)dt  =  A  =  limA(t). 

r— +00  T—*oo  ^  t—KX> 

Recall  that  T^  =  7.22176,  for  the  parameter  values  given  in  Section  6.3. 
This  value  is  entered  for  T  in  the  availability  function  in  Mathematica  and 
Aavgit),  is  calculated  for  t  =  500.,  750.,  and  1000.  These  results  are  shown  in 
Table  6.2  with  the  difference  from  the  limiting  availability,  A  —  0.734706,  noted 
in  the  "error"  column. 


Table  6.2  Limiting  Values  of  Aavg{t) 


t 

^avgit) 

Error 

500. 

0.737231 

0.002525 

750. 

0.736391 

0.001685 

1000. 

0.735970 

0.001264 
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As  expected,  the  error  decreases  as  t  increases.  At  f  =  1000.,  the  error  is 
approximately  10”^.  Given  the  approximations  employed  to  develop  the 
estimated  availability  function,  these  values  validate  the  performance  of  the 
estimated  function  as  time  increases  to  large  values. 


6.3. 1.3  Model  Recognition  of  Local  Maxima 

In  this  validation  test,  the  parameter  values  for  lifetime,  repair  and 
preventive  maintenance  densities  are  changed  to  the  values  shown  below; 

Weibull  lifetime  density;  a=  /3  =  2,  meanlife  =  12.5, 

characteristic  life  =  14.14 

Exponential  repair  time  density;  =  800.,  mean  =  0.00125 
Exponential  preventive  maintenance  density;  =  4000.,  mean  =  0.00025 

These  parameter  values  for  the  lifetime,  repair  time  and  preventive 
maintenance  time  densities  reflect  values  of  p  =  10000.  and  5  =  5.  Recall  that 
p  is  the  ratio  of  the  mean  lifetime  to  mean  repair  time  and  6  is  the  ratio  of  the 
mean  repair  time  to  mean  preventive  maintenance  time. 

These  parameter  values  are  entered  in  Mathematica  and  the  average 
availability  at  f  =  750.  is  evaluated  for  T  =  6.722,  7.222,  7.722,  and  8.222. 
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Mathematica  evaluates  these  average  availabilities  through  numerical 
integration.  Table  6.3  shows  the  results  obtained. 


Table  6.3  Local  Maxima  for  T^ 


T 

Aavgit  =  750.) 

6.722 

0.99992778 

T;.,  =  7.222 

0.99992801 

7.722 

0.99992793 

8.222 

0.99992762 

Due  to  the  large  value  of  p  used  in  this  validation,  availability  values  are 
very  high.  However,  even  with  such  a  large  ratio  of  mean  lifetime  to  mean 
repair  time,  the  estimated  availability  function  still  produced  a  local  maximum 
at  T*  =  7.222 . 

OC 


6.3.2  Availability  Model  Results 

For  the  following  plots,  the  Weibull  lifetime  density  function  is  used. 
Parameter  values  are  a=^  (scale  parameter)  and  p=2  (shape  parameter). 

The  parameter  values  provide  a  mean  component  lifetime  of  approximately 
12.5  and  a  characteristic  life  of  14.14.  The  repair  time  and  preventive 
maintenance  time  are  modeled  as  exponential  distributions.  The  mean  lifetime 
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to  mean  repair  time  ratio,  p,  is  set  to  a  value  of  p  =  2.  to  accentuate  the 
differences  in  the  availability  values.  The  value  of  the  mean  repair  time  to 
mean  preventive  maintenance  time,  6,  is  set  to  a  value  of  5  =  5.  Thus  the  mean 
repair  time  is  6.25  and  the  mean  preventive  maintenance  time  is  1.25 

In  the  following  series  of  plots,  the  first  plot  represents  the  estimated 
availability  function  for  0.1<f<15.  with  an  age  replacement  period  of 
T=7.22176.  This  T  value  is  the  optimal  replacement  period  for  an  infinite  time 
horizon.  Note  that  the  function  converges  rather  quickly  to  the  limiting 
availability  value  for  the  optimal  value,  T^.  This  value  of  A  is  approximately, 
>1=0.7347. 

The  subsequent  plots  show  the  estimated  availability  function  for  T=  1 , 

3,  5,  6.75,  12.5,  14.14,  28.28,  35,  42.4,  50.,  and  55.  as  a  function  of  time.  For 
reference  purposes,  the  T^=7.22176  plot  is  included  in  each  plot. 
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A(t) 


T=6.75 


Figure  6.5  Availability  Plots  for  T=6.75  and  1=1.22 


A(t) 


T=12.5 


t,  time 


Figure  6.6  Availability  Plots  for  T=12.5  and  T=7.22 
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A(t) 


T=14.14 


Figure  6.7  Availability  Plots  for  T=1 4. 14  and  T=7.22 


A(t) 


T=28.28 


t,  time 


Figure  6.8  Availability  Plots  for  T=28.28  and  T=7.22 
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Figure  6.12  Availability  Plots  for  T=55.  and  T=7.22 
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6.3.3  Availability  Model  Results  Discussion 


Figures  6.2,  6.3  and  6.4  show  an  anomaly  in  the  estimated  availability 
function,  A{t).  Note  the  portion  of  the  plot  that  goes  above  a  value  of  one  for 
A(t).  Since  the  availability  function,  A(t),  is  the  probability  that  the  component 
is  functioning  at  time  t,  the  maximum  A(t)  value  is  one.  This  error  is  small  and 
ranges  from  about  0.5  for  an  age  replacement  value  of  T=1 .  to  about  0.025  for 
an  age  replacement  value  of  T=5.  The  error  decreases  as  the  age 
replacement  value,  T,  increases.  Additionally,  the  error  occurs  later  in  the  time 
domain,  t,  as  the  age  replacement  value,  T,  increases.  This  small  anomaly 
occurs  due  to  the  estimated  Laplace  transform  inversion  technique  employed 
to  estimate  A(t).  Recall  that  the  methodology  employed  to  find  the  availability 
function,  A{t),  provided  an  exact  representation  of  the  Laplace  transform  of 
A{t).  The  only  estimation  technique  employed  is  the  use  of  the  Stehfest 
numerical  Laplace  transform  inversion  technique.  Thus,  this  estimated 
inversion  technique  caused  a  small  anomaly  to  occur  for  small  values  of  T. 

Figure  6.12  demonstrates  that  the  estimated  availability  function 
diverges  for  large  values  of  T  and  large  values  of  continuous  time,  f.  The 
limitation  of  the  model  is  T  <  55.  This  limitation  occurs  due  to  computational 
limitations  such  as  computer  round-off  error  and  floating  point  limitations. 

The  average  availability,  Aavgir),  for  a  finite  time  period  r,  may  be 
evaluated  graphically  from  Figures  6.1  through  6.1 1 .  Recall  that  Aavgir)  is 

Aavgir)  =  yfo  Aavgiu)du  . 
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For  a  finite  time  period,  r,  Aavgir)  may  be  graphically  evaluated  by  considering 
the  area  under  the  availability  function,  A{t),  curve  between  f  =  0  and  f  =  r, 
where  t  is  the  continuous  time  variable.  Furthermore,  if  a  comparison  of  the 
Aavgir)  for  two  different  functions  shown  on  the  same  plot  is  desired,  one  may 
simply  consider  only  the  intersection  of  the  areas  bounded  by  both  curves.  For 
example  consider  a  finite  time  period  of  1 0  time  units  and  refer  to  Figure  6.11. 
Note  that  for  the  time  period  from  t  =  0.  io  t  =  2.  that  the  curve  corresponding 
to  T^=  7.22  bounds  the  area  under  the  curve  for  T=50.  Since  the  T^=  7.22 
curve  bounds  a  larger  area  over  this  time  period,  T^=  7.22  has  a  larger 
average  availability  than  T=50.  Now  consider  the  time  period  from  f  =  2.  to 
t  =  10.  Note  that  in  this  case,  the  T=50.  curve  bounds  the  T^=  7.22  curve. 
Thus  for  this  time  interval,  T=50.  provides  a  larger  average  availability  value 
than  T^=  7.22  since  the  area  under  its  curve  is  greater.  In  fact  if  we  consider 
the  time  period  t  =  0Aot  =  10.,  it  is  graphically  apparent  from  Figure  6.1 1 ,  that 
the  T=50.  curve  bounds  more  area  than  the  T^=  7.22  curve.  Thus  the  age 
replacement  value  of  T=50.  provides  a  larger  value  of  Aavgir  =  10.)  than  the 
optimal  infinite  time  horizon  value,  T^=  7.22.  This  result  is  extended  further  in 
Chapter  7. 
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Chapter  7 


Conclusion  and  Accomplishments 


7.1  Availability  Model  Performance 

The  research  methodology  presented  produces  an  exact  representation 
of  the  Laplace  transform  of  the  availability  function,  A{t,  T).  This  availability 
function  represents  the  probability  that  a  component  functions  at  time  t  with  a 
preventive  maintenance  period  of  T  time  units.  T  represents  the  age  at  which 
the  component  is  preventively  maintained.  If  the  component  operates  T  time 
units,  it  is  preventively  maintained.  If  the  component  fails  before  T  time  units 
have  passed,  it  receives  a  failure  maintenance  action.  This  type  of 
maintenance  policy  is  commonly  known  as  an  age  replacement  preventive 
maintenance  policy.  However,  this  research  does  not  define  the  preventive 
maintenance  as  a  component  replacement.  This  research  requires  that  the 
preventive  maintenance  restore  the  component  to  as  good  as  new  condition. 
The  actual  action  taken  is  defined  by  the  specific  situation. 

The  derived  model  treats  the  age  replacement  period,  T,  as  a 
parameter.  Thus,  the  component  availability  at  a  given  time,  t,  is  dependent 
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on  the  age  replacement  period,  T.  The  methodology  presented  produced  the 
Laplace  transform  for  A(t,  T).  This  transform  could  not  be  inverted  exactly. 
Thus  an  inversion  estimation  technique  is  applied  to  obtain  an  estimate  of 
A{t,  T).  This  estimation  process  causes  some  anomalies  for  small  values  of  T. 
These  are  noted  in  Section  6.3.2.  Validation  tests  of  the  model  in  both  the 
exponential  case  as  well  as  the  Weibull-exponential  case  are  successful  and 
reviewed  in  Sections  7.1.1  and  7.1.2. 


7.1.1  Exponential  Validation  Case 

In  this  case,  the  component  lifetime  probability  distribution  as  well  as  the 
failure  repair  time  and  preventive  maintenance  service  time  are  modeled  as 
exponential  probability  distributions.  Each  distribution  is  considered  to  be 
unique  and  independent  of  each  other. 

This  case  is  considered  a  validation  case.  Recall  that  if  a  component's 
lifetime  is  considered  to  be  exponential,  then  no  benefit  will  result  from 
performing  preventive  maintenance.  Under  exponential  assumptions,  the 
component  does  not  exhibit  the  effects  of  aging.  No  matter  how  long  the 
component  has  been  in  operation,  the  probability  of  failure  in  the  next  instant 
of  time  remains  constant  (constant  failure  rate,  CFR).  Since  the  component 
does  not  age  under  these  assumptions,  no  benefit  results  from  purposely 
removing  an  operating  component  for  preventive  maintenance. 
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Three  cases  are  implemented  to  validate  the  availability  model  in  the 
exponential  case.  All  three  cases  support  the  validity  of  the  availability  model 
developed  and  implemented  in  this  research.  First,  the  model  is  tested  for 
performance  near  time  equal  to  zero  (f  =  0).  The  model  assumes  that  the 
component  is  working  at  f  =  0  and  the  availability  at  this  time  is  unity, 

A{t  =  0)  =  1.  Results  presented  in  Section  5.4  support  model  performance  in 
this  time  region. 

Second,  the  model  is  tested  at  very  large  values  of  time.  Values  of  the 
estimated  availability  function  at  these  large  values  of  time  should  converge  to 
the  limiting  availability  value.  Previous  research  results  provide  the  values  of 
the  limiting  availability,  A.  The  model  is  evaluated  at  large  values  of 
continuous  time,  t,  for  several  values  of  the  age  replacement  period,  T. 

Results  presented  in  Section  5.4  support  the  availability  model  performance  as 
the  continuous  time  variable,  t,  grow  large.  In  the  tests,  the  model  produces 
results  matching  the  limiting  availability  value  obtained  from  previous  research 
results. 

The  first  two  methods  of  validation  support  the  availability  model  validity 
as  t-*0  and  as  f->oo.  However,  model  validity  between  these  extremes  is  not 
supported.  The  third  test  provides  model  validation  for  times  between  these 
extremes.  This  test  takes  advantage  of  the  exponential  assumptions 
discussed  above.  If  the  component's  lifetime  is  assumed  to  be  exponential 
then  the  availability  model  should  produce  increasing  values  of  availability  as 
the  age  replacement  period,  T,  is  increased.  In  this  case,  maximum  average 
availability  as  well  as  limiting  availability  is  attained  by  performing  no 
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preventive  maintenance.  Thus,  the  age  replacement  period,  T,  should 
approach  infinity  (T-+00).  The  estimated  availability  function  produces  this 
result.  Table  5.1  shows  that  as  T  is  increased,  the  average  availability, 
^aug(f.T).  increases.  This  is  demonstrated  at  several  values  of  continuous 
time. 

Overall,  the  exponential  case  validated  the  availability  model  produced 
by  this  research.  The  model  exhibited  the  expected  behavior  at  extreme 
values  of  continuous  time,  t,  as  well  as  at  times  between  these  extremes. 


7.1.2  Weibull  -  Exponential  Case 

The  Weibull-exponential  case  assumes  that  the  component  lifetime 
probability  distribution  is  Weibull.  Furthermore,  the  shape  parameter  for  the 
Weibull  is  considered  to  be  greater  than  zero.  This  causes  the  component 
lifetime  distribution  to  have  an  increasing  failure  rate  (IFR).  Under  this 
assumption,  the  component  exhibits  the  effects  of  aging  since  the  probability 
that  the  component  fails  in  the  next  instant  of  time  is  now  dependent  on  the 
component's  age  or  in  effect  how  long  the  component  has  been  operated. 
Both  failure  repair  time  and  preventive  maintenance  time  probability 
distributions  are  modeled  as  unique  exponential  probability  distributions.  All 
probability  distributions  are  considered  to  be  independent  of  each  other. 
Three  validation  tests  are  performed  in  this  case.  The  model  is  tested  at  the 
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extreme  values  of  continuous  time  (f-»0,  t-»oD).  The  model  is  also  tested  for 
recognition  of  average  availability  local  maxima. 

In  the  Weibull-exponential  case,  results  validated  the  model 
performance  at  very  small  values  of  continuous  time,  f-»0.  Table  6.1  details  the 
results.  Values  of  the  estimated  availability  function,  A(t,  T)  are  taken  at 
t  =  0.000001  for  values  of  T  between  0.1  and  750.  In  all  cases,  a  value  of  unity 
is  obtained.  This  matches  the  expected  theoretical  value.  The  model  assumes 
that  the  component  is  functioning  at  time,  t  =  0.  By  definition,  the  availability  at 
f  =  0.  is  the  probability  the  component  is  functioning  at  f  =  0.  Thus,  the 
expected  availability  value  at  f  =  0.  is  one.  Note  that  the  estimated  availability 
function,  A(t,T)  for  the  Weibull-exponential  case  experiences  discontinuities  at 
t  —  0.  due  to  the  Laplace  inversion  estimation  and  the  Weibull  distribution 
transform  estimation  techniques  employed.  Therefore,  a  value  oft  =  0.000001 
is  used. 

Results  also  validate  the  model  performance  at  large  values  of 

continuous  time  (t-»oo).  For  this  case  the  optimal  replacement  period  for  the 

infinite  time  horizon,  T^,  is  calculated  per  previous  research  results.  The 

corresponding  limiting  value  of  availability,  A  is  also  calculated.  As  discussed 
earlier,  the  limiting  availability,  A  =  lim  A(f)  and  the  \im  Aavg{t)  are  the  same 

^--♦00  t— >oo 

when  the  limits  exist.  Thus,  the  values  of  the  estimated  Aavg{t,Y^)  should 
converge  to  A  for  large  values  of  t.  The  value  of  ylav5(f,T^)  is  calculated  for 
large  values  of  time  from  the  estimated  A{t,  T^).  Mathematica  is  used  to 
numerically  integrate  A(f,T^)  to  obtain  Ao„j(f,T^).  Table  6.2  details  the 
results  and  shows  an  error  of  1 0'^  at  f  =  1000. 
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The  third  validation  test  performed  for  the  Weibull-exponential  case 
tested  whether  the  model  could  recognize  an  average  availability  local 
maxima.  The  test  calculated  values  of  Aavg{t,  T)  for  values  of  T,  close  to  and 
including  the  value  of  T^  discussed  above.  In  this  validation  test,  limiting 
values  are  required,  thus,  a  large  value  of  time,  t  =  750.,  is  used.  Table  6.3 
details  the  results.  The  model  recognizes  T^  as  the  replacement  period  that 
maximizes  average  availability  for  very  large  values  of  time  {t-*oo). 

The  validation  tests  for  the  Weibull-exponential  case  demonstrate 
model  performance.  Results  validate  model  performance  at  small  values  of 
time  (<->0)  and  at  large  values  of  time  (t-^oo).  Also,  the  model  is  shown  to 
recognize  T^  as  the  age  replacement  period  that  maximizes  average 
availability  for  very  large  values  of  time  (t-*oo). 


7.2  Maximizing  Average  Availability  for  Finite  Component  Economic 
Lifespan 

The  results  in  Section  6.3.2  show  that  for  finite  time  intervals,  the  optimal 
infinite  time  horizon  replacement  interval,  T^,  does  not  maximize  average 
availability.  This  result  is  shown  for  example  in  Figure  6.1 1  and  is  discussed  in 
Section  6.3.2.  The  conclusion  is  that  the  optimum  infinite  time  horizon 
preventive  maintenance  policy,  T^,  does  not  maximize  average  availability  for 
some  finite  time  periods. 
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In  actual  operation,  repairable  components  and  systems  are  designed  to 
be  used  over  a  finite  time  interval.  These  lifecycle  design  constraints  are  very 
important  in  determining  the  component  and/or  system  maintenance 
requirements.  This  research  concludes  that  using  the  age  replacement  period 
based  upon  an  infinitely  long  operational  interval  does  not  maximize 
component  and/or  system  availability  for  a  finite  operating  design  lifespan. 
This  result  is  important  to  system  and/or  component  lifecycle  design 
considerations. 


7.3  Summary  of  Significant  Research  Accomplishments 

This  research  presents  a  unique  and  original  approach  to  developing 
the  component  availability  function  for  the  age  replacement  preventive 
maintenance  policy.  The  approach  explicitly  considered  the  effects  of  unique 
failure  repair  time  and  preventive  maintenance  action  time  probability 
distributions.  The  availability  model  included  derivation  of  the  partial  Laplace 
transform  for  cases  where  the  time  variable  is  truncated  (i.e.  does  not 
approach  infinity).  Normal  Laplace  transforms  require  the  time  variable  to 
have  a  range  of  [O.oo].  The  end  result  is  an  estimate  of  component  availability 
as  a  function  of  time,  t,  as  well  as  the  age  replacement  period,  T.  Thus,  the 
component  availability  function  is  defined  as  A(i,T). 

The  availability  model  derived  in  this  research  produces  an  exact 
representation  of  the  Laplace  transform  the  component  availability  function. 
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This  is  defined  as  ^’(a,T).  Even  in  the  exponential  case,  could  not  be 

inverted  exactly  to  find  i4(t,T).  A  numerical  Laplace  transform  inversion 
estimation  technique,  the  Stehfest  method,  is  employed.  The  Stehfest  method 
and  other  inversion  estimation  techniques  have  been  used  to  find  estimates  of 
the  time  domain  function  at  specific  values  of  time.  Thus,  past  research 
applied  these  techniques  for  numerical  results.  However,  this  research 
extended  the  use  of  the  Stehfest  method  to  estimating  the  function  A{t,  T) 
from  A"{s,  T).  The  Wolfram  Research,  Inc.  software  package,  Mathematica, 
allowed  the  Stehfest  series  estimation  technique  to  be  applied  in  general  and 
thus  produce  an  estimate  of  the  availability  function,  A{t,  T). 

In  the  Weibull-exponential  case,  the  partial  Laplace  transform  must  be 
evaluated  for  the  Weibull  distribution  and  the  Weibull  survivor  function.  Recall 
that  in  the  partial  Laplace  transform,  the  time  variable  is  truncated  and  thus 
ranges  over  [0,  r],  where  r  is  finite.  This  research  developed  an  exact  infinite 
series  representation  for  the  partial  Laplace  transform  of  the  Weibull 
distribution  and  survivor  function.  This  development  also  yielded  the  normal 
Laplace  transform  for  the  Weibull  distribution  and  survivor  function.  In  the 
normal  Laplace  transform  the  time  variable  ranges  over  [0,  oo].  The  latter 
result  is  particularly  useful  in  many  applied  probability  applications  such  as 
reliability  studies. 

Finally,  this  research  demonstrated  that  average  availability  for  a  finite 
time  periods,  r,  is  not  maximized  by  the  optimal  infinite  time  horizon 
replacement  period,  T*^,  at  all  values  of  r.  This  result  is  important  to  lifecycle 
design  considerations  for  maintaining  and  servicing  components  and  systems. 
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Chapter  8 


Further  Research 


8.1  Finite  Economic  Life  Model 


In  the  availability  model  developed  and  presented  in  this  research  the 
component's  economic  life  is  not  bounded.  This  model  is  useful  for  estimating 
the  effects  of  the  infinite  time  horizon  age  replacement  period  (preventive 
maintenance  policy),  T^,  on  a  finite  economic  component  life.  Definitive 
values  of  the  age  replacement  period  value,  T,  require  further  research  and 
development  of  the  availability  model. 

Let  T  be  defined  as  the  economic  life  of  the  component.  Equation  (4.8)  is 
modified  to 


A{t) 


RL{i)+fo  m{u)  du,  for  f  <  T, 


Ri{t—u)  m{u)  du, 


for  T<  f  <  T. 


(8.1) 


Note  that  the  continuous  time  variable,  t,  is  no\A^  bounded  by  the  economic 
component  life,  r.  As  in  Chapter  4,  the  development  proceeds  with  the 
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construction  of  the  Laplace  transform  of  the  availability  function,  A(f),  stated  in 
equation  (8.1 ).  Recall  that  the  definition  of  the  Laplace  transform  for  the 
availability  function  is 


A*(s)  =  fJ°e-^^A(t)dt. 

Applying  (8.1)  to  this  relationship,  the  following  is  obtained: 

A*(s)  =  -h  Jq  u)m(u)du^dt^ 

+  ^  Jj  e~^ *'  ^ Riit  —  u)  m(u)  dv^  d^ 

—  lo  ^Lii)dt  +  (^jQRi{t  -  u)m(u)  dujdt 

+  Jj  Ri{t  —  u)m(u)  dxi^dt 

—  Jq  (l-Fiit))  dt  +  JJJq  Ri(t  -  u)  m{u)  du  dt 

+  fj  Ri^t  —  fn(u)  du  dt  . 


116 


=  fj  dt  -  e-^^FUt)  dt  +  fjfj  -  n)  m(u)  dt  du 

+  Jo  Jj*'"  Riit  -  u)m{u)dtdu  +  JJ  Jle~^*- Ri{t- u)m{u)dtdu 

=  -  FI{3,  T)  +  fj  m(u)  fj  Riit  -  u)  dt  du 

+  fo  Jt*'^  Riii  —  u)  dtdu  +  JJ  m{u)  fJ  Ri{t  —  u)  dtdu 

=  -  ^1(3,  T)  +  fJ  m(u)  Ri(y)  dy  du 

+  fJ  m(u)  +  RL(y)  dydu  +  fJ  m(u)  /(["“  +  Rziy)  dydu 

=  -  Fl{s,  T)  +  fJ  e— m(ti)  /J-“  Ri(y)  dy  du 

+  Jo  e"^“"z(u)  /7_u  e“®2/  dydu  +  fJ  e“^“m(u)  /J'““  e”"*'  R^(y)  dydu 

=  ^  -  ^2(«.  T)  +  fJ  e— m(tx)  fJ  e-^y  R^iy)  dy  du 

+  /f  e'^  “m(u)  /J’"  “  e”^  "  i?i  (y)  dy  rfu 
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=  (‘^  -  T))  (l  +  lo  e-’-miu)  du) 

+  It  e~^“m(u)  /J"""  Ri{y)  dydu 

+  /t  e“"S'  (1  -  Fiiy))  dydu 

=  {^-Fl{s,  T))  (l+m*(s,  T)) 

+  fje~^^m(u)  -  Fl(s,  T-u)jdu 

={^-Flis,J))(l+m*{s,T)) 

+  It  “  It  F2(s,  t  —  u)du  (8.2) 
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The  second  term  of  the  sum  shown  in  equation  (8.2)  may  be  restated  as 
follows: 


It  (-•■■"f-'-)  du  =  m(u)  du-^  Jj  m{u)  du 

=  ;  (rn*(s,  r)  -  m’(s,  T)j  -  ^  (^fjm(u)  du  -  fjm(u)  du^ 

=  i  (rn*(s,r)  -  m*(s,T)J  -  -  M(T)j 

Substituting  this  result  into  equation  (8.2),  we  obtain 

-Fl{s,  T))  ^l  +  m*(s,T))  +  i(Tn*(s,r) -Tn*(s,T)) 

-  ^  (m(t)  -  M(T))  -  /;  e-“m(ti)  F£(s,  T-u)du  (8.3) 

Equation  (8.3)  provides  the  Laplace  transform  of  the  availability  function  for  a 
finite  economic  life  span,  r.  The  renewal  density,  m(it),  is  defined  and 
discussed  in  Section  4.1.1.  This  density  function  may  be  found  through  the 
density  function  for  the  time  between  renewal,  h{t).  M{u)  is  the  renewal 
function  and  is  defined  as 

M(u)  =  /;  m{x)  dx. 

Fl{s,  T)  is  the  partial  Laplace  transform  of  the  component's  lifetime  probability 
distribution. 
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Additional  research  of  this  model  is  required.  Methods  to  evaluate  this 
model  need  to  be  determined.  This  includes  numerical  approaches  to  evaluate 
the  availability  function  as  well  as  optimization  techniques  to  find  an  optimal 
value  of  T  that  maximizes  average  availability  over  the  component's  finite 
economic  life  span.  One  specific  complication  is  the  requirement  for  the 
renewal  density  and  renewal  function.  Numerical  methods  to  accurately  find  or 
estimate  these  functions  must  be  found. 


8.2  Service  Time  Probability  Distributions 

The  specific  cases  examined  in  this  research  assume  that  the  repair 
time  and  the  preventive  maintenance  service  time  probability  distributions  are 
exponential.  Section  2.3.3  provides  a  reference  for  this  assumption  [88]. 
However,  other  distributions  may  more  accurately  model  the  randomness  of 
service  times  since  the  exponential  distribution  is  "memoryless".  Under 
exponential  assumption,  if  a  component  is  serviced  for  z  time  units,  then  the 
probability  that  the  unit  is  repaired  in  the  next  instant  of  time  remains  constant 
for  all  values  of  2.  Essentially,  this  means  that  if  you  have  worked  on  repairing 
(or  preventively  maintaining)  a  component  for  1  second,  1  minute,  1  hour,  1 
day,  1  month  or  1  year  then  the  probability  that  you  fix  the  component  in  the 
next  instant  of  time  remains  constant  no  matter  how  long  you  have  worked.  In 
most  practical  applications,  the  probability  the  component  is  fixed  increases  as 
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the  repair  time  increases.  Thus  a  more  practical  probability  distribution  will  be 
IFR  (increasing  failure  rate). 

Two  well  known  distributions  may  be  IFR.  The  Weibull  distribution  is  IFR 
for  values  of  the  shape  parameter  greater  than  1  (y0  >  1 ).  The  log-normal 
distribution  may  be  both  IFR  and  DFR  (decreasing  failure  rate)  [95].  For 
example,  for  service  times  less  than  some  threshold  value,  the  log-normal 
distribution  may  be  DFR  and  for  service  times  above  this  value  it  may  be  IFR. 
However,  the  use  of  either  distribution  presents  complications.  The  model 
developed  in  this  research  requires  the  Laplace  transform  of  the  repair  time 
and  preventive  maintenance  service  time  distributions.  This  research  presents 
the  Laplace  transform  for  the  Weibull  distribution.  However,  the  Laplace 
transform  for  the  log-normal  distribution  needs  to  be  developed.  The  Weibull 
Laplace  transform  is  an  infinite  series  and  it  is  highly  probable  that  the  log¬ 
normal  may  also  be  an  infinite  series. 

In  either  case,  the  infinite  series  must  be  truncated  to  estimate  the 
required  Laplace  transform.  Considering  equation  (4.16)  this  may  require  the 
multiplication  of  two  expressions  containing  many  terms.  For  example, 
consider  a  case  where  we  have  Weibull  component  lifetimes  with  unique 
Weibull  repair  and  preventive  maintenance  times.  In  the  Weibull-exponential 
case,  32  terms  were  used  for  the  partial  Weibull  Laplace  transform.  If  we  use 
32  terms  in  this  new  case,  we  are  faced  with  multiplying  a  32  term  expression 
by  another  32  term  expression.  This  results  from  multiplying  fl{s,T)  and  g*{s). 
This  represents  the  first  term  in  equation  (4.16).  Obviously,  the  complexity  has 
increased.  Methods  to  reduce  the  complexity  need  to  be  found. 
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APPENDIX  1 

Numerical  Laplace  Transform  Inversion  Test 

The  purpose  of  this  appendix  is  to  demonstrate  that  Mathematica  allows 
the  estimate  of  a  time  domain  function  through  the  use  of  the  numerical 
Laplace  transform  inversion  package,  niapinv.m  .  This  numerical  Laplace 
transform  inversion  package  is  available  from  the  Wolfram  Research,  Inc. 
world  wide  web  (WWW)  site,  http;//www. wri.com. 

As  explained  in  the  main  text,  the  numerical  inversion  techniques  were 
developed  to  obtain  estimates  of  the  time  domain  function  at  specific  numerical 
values  of  time.  These  numerical  Laplace  transform  inversion  techniques  were 
developed  to  be  applied  to  Laplace  transforms  that  cannot  be  inverted  exactly 
through  conventional  mathematical  techniques.  The  Stehfest  method 
estimates  this  value  through  the  application  of  a  mathematical  formula 
containing  two  levels  of  summations.  By  leaving  the  time  specification  as  a 
parameter,  t.  in  this  case,  Mathematica  produces  an  estimate  of  the  time 
domain  function  (as  a  function  of  t)  by  evaluating  the  series  contained  in  the 
Stehfest  formula  explicitly.  This  methodology  results  in  an  estimate  for  the 
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time  domain  function  containing  over  100  terms  for  N=6.  For  this  reason,  the 
time  domain  function  is  not  displayed. 

The  first  Mathematica  statement  defines  the  test  function  A(s).  In  this 
case  it  is  equation  (5.3).  Mathematica  input  statements  are  represented  by 
boldface  type.  Mathematica  output  resulting  from  the  input  statements  is 
represented  by  italicized  type. 

A[sJ=  ((P+s)  (R+s))/(((P+s)  (R+s)  (L+s))-(L  R  (P+s))+ 

((L  R  (P+s))-{P  (L+s)  (R+s)))  Exp[-T  (L+s)]) 

((P  +  s)  (R s))  / 

(-(L  R(P  +  s))  •i‘(L  +  s)(P  +  s)  (R+s)  + 

LR(P  +  s)-P(L  +  s)(R+s) 

(L+s)T 

E 

The  following  statement  loads  the  numerical  Laplace  transform 
inversion  package  available  from  Wolfram  Research,  Inc.  In  this  case  the 
Mathematica  file  containing  the  package  is  'niapinv.m'  and  is  located  in  the 
'wnmath22’  subdirectory  on  drive  c. 

«c:\wnmath22\nlapinv.m 

The  following  statement  shows  the  Mathematica  syntax  for  using  the 
Stehfest  method  for  numerical  Laplace  transform  inversion. 
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?NLInvSteh 


NUnvSteh[expr,  s,  t,  n] 

where  expr  =  the  Laplace  transform  expression  to  be  inverted 
s  =  the  Laplace  transform  parameter 
t  =  time 

n  =  number  of  terms  in  series  (Hi  must  be  an  even  number  Hi), 
typically  6  <n  <20. 


The  function,  lnvA1  (t),  is  defined  as  the  Stehfest  numerical  Laplace 
transform  inversion  of  the  Laplace  transform  function,  A(s).  Recall  that 
lnvA1  (t)  represents  an  estimate  of  the  time  domain  function  obtained  by 
numerically  inverting  A(s). 


lnvA1  [t J=NLInvSteh[A[s],s,t,6]; 


The  variables  L,  R  and  P  define  the  exponential  probability  density 
parameters  for  the  lifetime,  repair  time  and  preventive  maintenance  time 
densities  respectively.  These  parameters  were  not  numerically  defined  when 
A(s)  was  entered.  The  parameters  are  numerically  defined  below.  The  mean 
lifetime  to  mean  failure  repair  time  ratio  is  1000  to  1  and  the  mean  failure 
repair  time  to  mean  preventive  repair  time  ratio  is  2.5  to  1 . 

L=1. 

R=1000. 

P=2500. 

1. 

1000. 

2500. 
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The  age  replacement  parameter,  T,  is  set  to  0.6  for  numerical 
evaluation.  Recall  that  T=0.6  means  that  a  component  is  preventively 
maintained  when  it  attains  an  age  of  0.6  time  units.  Complete  renewal  is 
assumed  for  either  failure  or  preventive  maintenance  actions  thus  the 
component's  age  is  reset  to  zero  after  either  action  is  completed. 


T=0.6 

0.6 

Evaluating  lnvA1(t)  at  t=5.,  we  obtain  a  time  domain  function  value  of 
2.19671. 

lnvA1[5.] 

2. 19671 

The  following  statement  takes  the  original  Laplace  transform  function, 
A(s),  and  directly  inverts  the  function  numerically  using  the  Stehfest  method  at 
t=5.  The  value  obtained  is  exactly  the  same. 

NLInvSteh[A[s],s,5.,6] 

2.19671 

The  age  replacement  parameter,  T,  is  reset  to  a  variable  through  the 
Mathematica  "Clear"  statement. 
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Clearrn 


The  following  two  plots  demonstrate  that  lnvA1(t)  is  equivalent  to  the 
direct  Stehfest  numerical  inversion  of  A(s)  at  discrete  points  in  time.  In  the  first 
plot,  the  estimate  time  domain  function,  lnvA1  (t),  is  plotted  for  values  of  t 
between  0.1  and  15.  and  T  between  0.05  and  10.  In  the  second  plot,  the 
ranges  of  the  variables  t  and  T  are  exactly  the  same.  However,  note  that  the 
function  evaluated  is  the  direct  Stehfest  numerical  inversion  of  the  original 
Laplace  transform  function,  A(s). 

Plot3D[lnvA1  ,1 5.},{T,0.05,1 0.}] 


15 


Fig.  A1.1  Plot  of  InvAI(t),  0.1  <  t  <  15.,  0.05  <  T  <  10. 
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Plot3D[NLInvSteh[A[s],s,t,6],{t.0.1,15.},{T,0.05,10.}] 


15 

Fig.  A1.2  Plot  of  Direct  Stehfest  Inversion,  0.1  <  t  <  15.,  0.05  <  T  <  10. 

The  following  output  further  demonstrates  that  InvAI(t)  is  equivalent  to 
the  direct  Stehfest  numerical  inversion  of  the  original  Laplace  transform 
equation,  A(s).  The  'Table"  function  provides  specific  function  values  for  both 
InvAI(t)  and  NLInvSteh[A[s],  s,t,6].  For  the  ranges  of  t  and  T  evaluated,  the 
results  match  exactly. 


Table[lnvA1[t],{t,0.1,5.1,0.5>,{T,0.1,5.1,0.5}] 

{{1.42798,  1.0046,  0.999161,  0.999004,  0.999001,  0.999001,  0.999001, 
0.999001,  0.999001,  0.999001,  0.999001), 

{5.03848,  1.28497,  0.989437,  0.981246,  0.99531,  0.999234,  0.999514, 
0.999277,  0.999116,  0.999045,  0.999017), 

(7.19112,  1.66612,  1.18129,  1.03008,  0.993293,  0.99171,  0.995391, 
0.997742,  0.998701,  0.998994,  0.99905), 

{8.45334,  1.87848,  1.30932,  1.11227,  1.03173,  1.00353,  0.99663, 
0.996403,  0.997421,  0.998233,  0.998681), 
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{9.20271,2.00367,  1.38155,  1.16618,  1.06857,  1.02401,  1.00579, 
0.999572,  0.998077,  0.998084,  0.998398}, 

{9.65807,  2.07965,  1.42417,  1.19837,  1.09374,  1.04148,  1.01621, 
1.00496,  1.00049,  0.998994,  0.998649), 

{9.94162,  2.12695,  1.45032,  1.21786,  1.1098,  1.05396,  1.02486, 
1.0103,  1.00346,  1.00049,  0.99933), 

{10.1222,  2.15708,  1.46684,  1.22999,  1.11999,  1.06242,  1.03129, 
1.01472,  1.00623,  1.00209,  1.00019), 

{10.2394,  2.17664,  1.47752,  1.23773,  1.12653,  1.06807,  1.03588, 
1.01812,  1.00855,  1.00356,  1.00106), 

{10.3168,  2.18954,  1.48454,  1.24277,  1.13079,  1.07185,  1.03909, 
1.02065,  1.01038,  1.0048,  1.00185), 

{10.3685,  2.19817,  1.48922,  1.24609,  1.13361,  1.07441,  1.04134, 
1.0225,  1.01179,  1.0058,  1.00252}) 


Table[NUnvSteh[A[s],s,t,6],{t,0.1,5.1,0.5},{T,0.1,5.1,0.5}] 

{{1.42798,  1.0046,  0.999161,  0.999004,  0.999001,  0.999001,  0.999001, 
0.999001,  0.999001,  0.999001,  0.999001), 

{5.03848,  1.28497,  0.989437,  0.981246,  0.99531,  0.999234,  0.999514, 
0.999277,  0.999116,  0.999045,  0.999017), 

{7.19112,  1.66612,  1.18129,  1.03008,  0.993293,  0.99171,  0.995391, 
0.997742,  0.998701,  0.998994,  0.99905), 

{8.45334,  1.87848,  1.30932,  1.11227,  1.03173,  1.00353,  0.99663, 
0.996403,  0.997421,  0.998233,  0.998681), 

{9.20271,  2.00367,  1.38155,  1.16618,  1.06857,  1.02401,  1.00579, 
0.999572,  0.998077,  0.998084,  0.998398), 

{9.65807,  2.07965,  1.42417,  1.19837,  1.09374,  1.04148,  1.01621, 
1.00496,  1.00049,  0.998994,  0.998649), 

{9.94162,  2.12695,  1.45032,  1.21786,  1.1098,  1.05396,  1.02486, 
1.0103,  1.00346,  1.00049,  0.99933), 

{10.1222,  2.15708,  1.46684,  1.22999,  1.11999,  1.06242,  1.03129, 
1.01472,  1.00623,  1.00209,  1.00019), 

{10.2394,  2.17664,  1.47752,  1.23773,  1.12653,  1.06807,  1.03588, 
1.01812,  1.00855,  1.00356,  1.00106), 

{10.3168,  2.18954,  1.48454,  1.24277,  1.13079,  1.07185,  1.03909, 
1.02065,  1.01038,  1.0048,  1.00185), 

{10.3685,  2.19817,  1.48922,  1.24609,  1.13361,  1.07441,  1.04134, 

1.0225,  1.01179,  1.0058,  1.00252)} 
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APPENDIX  2 


Mathematica  Statements  for  the  Exponential  Case 


The  purpose  of  Appendix  2  is  to  provide  an  example  of  the  Mathematica 
input  statements  used  for  analysis  of  the  exponential  case  referenced  in 
Section  5.3.1 .  Note  that  Mathematica  input  statements  appear  in  boldface 
type  while  Mathematica  output,  except  graphic  output,  appears  in  italicized 
type. 

The  following  statement  loads  the  "DiracDelta"  package  in  Mathematica. 
This  package  defines  the  diracdelta  function  also  known  as  the  Unitstep 
function  for  use  in  Mathematica. 

«Calculus'DiracDelta' 

The  following  statement  defines  the  equation  (5.3)  which  is  also  the  first 
part  of  equation  (5.1 )  referenced  in  Section  5.3. 
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A[sJ=  ((P+s)  (R+s))/(((P+s)  (R+s)  (L+s)HL  R  (P+s))+((L  R  (P+s))-(P  (L+s) 
(R-:  s)))  Exp[-T  (L+s)]) 

((P  +  s)  (R  +  s))  / 

(-(L  R(P  +  s))  +(L+s)(P  +  s)  (R+s)  + 

LR(P  +  s)-P(L+s)  (R  +  s) 

- 

(L^s)T 

E 


The  following  statement  loads  the  numerical  Laplace  transform 
inversion  package,  "niapinv.m".  Note  that  the  package,  "nlapinv.m",  is  located 
in  the  "wnmath22"  subdirectory  on  drive  c. 

«c:\wnmath22\nlapinv.m 


This  statement  provides  the  correct  syntax  and  use  of  the  Stehfest 
numerical  inversion  package. 


?NLInvSteh 

NLInvStehfexpr,  s,  t,  n] 

where  expr  =  the  Laplace  transform  expression  to 
be  inverted 

s  =  the  Laplace  transform  parameter 
t  =  time 

n  =  number  of  terms  in  series  (HI  must  be 
an  even  number  HI), 
typically  6  <n  <20. 
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The  following  statement  defines  the  function,  lnvA1(t).  This  function  is  a 
numerical  inversion  of  the  function  A(s)  defined  above.  Note  that  the  Stehfest 
method  is  used  for  the  numerical  inversion.  The  resulting  function,  lnvA1(t),  is 
a  function  of  time,  t  and  the  age  replacement  parameter,  T. 

InvAl  [t J=NLInvSteh[A[s],s,t,6]; 


To  obtain  a  functional  estimate  of  the  availability  function,  equation  (5.2) 
is  applied  as  noted  in  Section  5.3.1 .  Note  in  this  statement  as  well  as  the 
previous  statement,  the  actual  output  of  the  functions  are  suppressed  by  the 
semicolon  appearing  at  the  end  of  the  statement.  These  functions  are  created 
by  applying  the  Stehfest  series  expansion  noted  in  Section  5.3.1 .  and  become 
very  long  due  to  nested  summations  within  the  Stehfest  method. 

TAvall[tJ=lnvA1[t]-(Exp[-T  L]  UnitStep[t-T]  InvAl  [t-T]); 

The  following  statements  set  the  exponentially  distributed  lifetime, 
failure  repair  time,  and  the  preventive  maintenance  time  parameters  to  1.,  10. 
and  25.  The  ratio  of  mean  lifetime  to  mean  repair  time  is  10  and  the  ratio  of 
mean  repair  time  to  mean  preventive  maintenance  time  is  2.5. 
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L=1. 

R=10. 

P=25. 


1. 

10. 

25. 


The  following  statement  produces  a  three  dimensional  plot  of  the 
availability  function  for  the  exponential  case  for  time,  t,  range  of  [0.1  20.]  and 
age  replacement,  T,  range  of  [0.05,  15.] 

Plot3D[TAvail[t],{t,0.1,20.},{T,0.05,15.}] 


Fig.  A2.1  A(t,T),  0.1  <  t  <  20..  0.05  <  T  <  15. 


The  following  two  dimensional  plots  represent  "slices"  of  the  three 
dimension  plot  shown  above.  The  variable  representing  continuous  time,  t,  is 
held  constant  while  the  age  replacement  period  is  allowed  to  vary.  The 


144 


response  is  the  availability  at  a  given  point  in  time  and  is  plotted  against  the 
age  replacement  variable.  Thus  the  graphs  shown  below  are  "slices"  of  Fig. 
A2.1  along  the  y-axis  (representing  the  age  replacement  period.  T)  for  constant 
points  in  time. 

Plot[TAvail[0.75], {1,0.001,0.749}] 


Fig.  A2.2  A(t,T).  t=0.75.  0.001  <  T  <  0.75 

Plot[TAvail[1.],{T,0.01,2.5}] 


Fig.  A2.3  A(t,T).  t=1.,  0.01  <  T  <  2.5 
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Plot[TAvail[2.],{T,0.1,3.}] 


Fig.  A2.4  A(t,T).  t=2..  0.1  <  T  <  3. 


Plot[TAvaiI[5.],{T,0.1,6.}] 


Fig.  A2.5  A(tJ).  t=5.,  0.1  <  T  <  6. 


Plot[TAvall  [1 0.],{T,0.1 ,7.}] 


Fig.  A2.6 


A(t,T).  t=10.,  0.1  <T< 


Plot[TAvail[20.],{T,0.1,10.}] 


Fig.  A2.7  A(t,T),  t=20.,  0.1  <  T  <  10. 
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Plot[TAvail[50.],{T,0.1,25.}] 


0.909 

0.90875 

0.9085 

0.90825 


0.90775 


15 


20 


Fig.  A2.8  A(t.T),  t=50.,  0.1  <  T  <  25. 
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